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Abstract 

Elastic interactions arising from a difference of lattice spacing between two co- 
herent phases can have a strong influence on the phase separation (coarsening) 
of alloys. If the elastic moduli are different in the two phases, the elastic inter- 
actions may accelerate, slow down or even stop the phase separation process. If 
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the material is elastically anisotropic, the precipitates can be shaped like plates 
or needles instead of spheres and can form regular precipitate superlattices. 
Tensions or compressions applied externally to the specimen may have a strong 
effect on the shapes and arrangement of the precipitates. In this paper, we re- 
view the main theoretical approaches that have been used to model these effects 
and we relate them to experimental observations. The theoretical approaches 
considered are (i) 'macroscopic' models treating the two phases as elastic media 
separated by a sharp interface (ii) 'mesoscopic' models in which the concentra- 
tion varies continuously across the interface (iii) 'microscopic' models which use 
the positions of individual atoms. 

Keywords: kinetics of phase separation, quenched alloys, elastic interactions, 
sharp interface model, diffuse interface models, atomic lattice models 
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1 Introduction 



1.1 The scope of this paper 

The separation of an alloy into two phases is a technologically important phe- 
nomenon whose theory owes much to the work of John Cahn. When the two 
phases first form they are very finely mixed, with large amount of interface be- 
tween the two phases. Atoms then diffuse within the material so as to reduce 
the amount of interface, making the mixture of phases more and more coarse 
as time proceeds. The driving force for this coarsening process comes from the 
free energy of the interface. For many alloys, where the sizes of the two kinds 
of atoms composing the alloy are not too different, the atoms manage to fit 
themselves on to a common lattice, but they can achieve this only at the cost 
of significant distortion of the lattice. At first, the energy arising from this 
distortion is relatively unimportant, but as the coarsening proceeds the intcrfa- 
cial energy decreases and eventually the elastic interactions, which are of long 
range and can also be strongly anisotropic, may predominate. This change in 
the dominant driving force can dramatically affect both the rate of coarsening 
and the structure of the domains themselves in the late stages of the coarsening 
process. 

In the present article we review some theoretical approaches that have been 
used to model the effects of elastic interactions on coarsening in alloys. We shall 
confine ourselves to alloys which are coherent (meaning that the lattice structure 
is only distorted, not disrupted, by the misfit between the two types of atom), 
and for that reason we shall not consider any effects involving dislocations, such 
as plasticity. Moreover, we shall consider only first-order phase transformations, 
in which the concentrations are different in the two phases, so that diffusive 
transfer of material from one place to another is an essential part of the process; 
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martensitic transformations, for example, will not be considered. 

An idea of the relative importance of the elastic and surface energy effects 
can be obtained by comparing the surface energy of a spherical precipitate of 
radius R, which is AiraR 2 where a is the interfacial energy per unit area (surface 
tension), with the elastic energy due to such a precipitate, which is of order 
i(47ri? 3 /3)Gq 2 where G is a typical elastic modulus, e.g. the shear modulus, 
and q is a typical strain due to the misfit. The two energies are equal when 
R ~ i?o = Qa/Gq 2 . At first, the precipitates are much smaller than this size, so 
that the interfacial tension predominates; but later on in the coarsening process 
the precipitates will be of this size or larger, so that elastic effects become 
important. As an example, the value of Ro for M3AI precipitates in nickel is 
of order 10 -7 m, and it has been shown [|43[ by experiments on a sequence of 
alloys with different values for Rq , that some of the effects due to elastic misfit 
are proportional to R/Ro at small R/Rq. 

We shall describe three different types of model which have been used in 
the theory of these elastic effects. In the first, the two phases are treated as 
continuous materials obeying the macroscopic laws of linear elasticity theory 
[], and the interface between them is treated as a geometrical surface (i.e. its 
thickness is ignored). This was the model used by J. D. Eshelby in his pioneering 



work 205 , [ 202 1, |196| . We shall call it the sharp interface model. There are 
two versions of this model : a 'static' version in which the energies of different 
arrangements of the phases are compared without detailed consideration of how 
the system can get from one arrangement to the other, and a 'dynamic' version 
in which the mechanism for such changes, namely diffusion, is also included in 
the model. The static version will be considered in chapter || and the dynamic 
version in chapter ^. 



1 Grinfeld [L3S| is alone in taking seriously the possibility of nonlinear effects 
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A second type of model, introduced by John Cahn in 1961, takes the struc- 
ture of the interface into account by using the concentration of one of the alloy 
components as a field variable. Such models are sometimes described as 'meso- 
scopic'. The concentration field is approximately constant (but at different 
values) in the two phases, and it varies continouously across the interface. Its 
time variation is given by a deterministic differential equation, analogous to the 



Cahn-Hilliard equation [ 195 1 . We shall call this the diffuse interface model and 
describe it in chapter |]. 

In the third type of model a completely microscopic description is used. 
The atoms are no longer represented by a continuum; instead we follow them 
individually and model their diffusion by random jumps. Although it gives no 
analytic results, this method has the merit that the physical assumptions going 
into it are very simple. It also takes fluctuations into account and therefore 
includes the possibility of nucleation in a natural way. We shall call it the 
atomic lattice model and describe it in chapter [5]. 

In cases where the lattice parameters are approximately the same in the 
two phases, so that the elastic misfit is unimportant, all three of these methods 
have been found to give results in reasonable agreement with experiments. For 



reviews, see Gunton|147|, Furukawa [ 144 1 , 13£], Binder [p3f. However, when the 
misfit is important, the theoretical situation becomes more complicated, and it 
may be that a combination of all the three methods, or more, is necessary for a 
proper understanding. 

Some effects on morphology and kinetics typically observed in experiments 
are summarized in the next section. Many of these effects are actually repro- 
duced by the various theoretical approaches. These effects range from a change 
(typically a slowing down) in the rate of precipitate coarsening to a transforma- 
tion from nearly spherical to cube-like or plate-like precipitates aligned along 
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the elastically soft directions, that is, [100] and equivalent ones in the case of 
Ni-base alloys. One also observes a direct effect of external stresses on the 
coarsening behaviour of alloys, which can lead to the formation of cylinders or 
parallel plates with wavy interfaces 

1.2 Experimentally observed effects due to elastic misfit 
interactions 

The aim of the present chapter is not to give an exhaustive review of the effects 
of elastic misfit strain on the phase separation process but rather to show some 
typical examples. Because of the huge practical importance of nickel-base su- 
peralloys for use at high temperatures, e.g. in turbine blades, these alloys have 
attracted much interest in recent experimental studies. Most of the following 
examples will focus on this type of alloy system where the nickel-rich matrix 
(7-phase) has a face-centered cubic structure and contains 7'-precipitates with 
Li2 structure of the type N13X, where X may be Al, Ti, Si, etc. 

1.2.1 Morphological effects 

Depending on the composition of the alloy, the precipitates may be either round 
as, e.g. in Ni-Al-Si |5(| or cube-like as in, for example, Ni-Al |36|. The cube- 
like shape is, in fact, very frequent in this type of system, appearing also in 
ternary alloys like Ni-Al-Mo |o| H, || H, 0, Ni-Al-Ti jiioj] , Ni-Al-Cr 



or Ni-Al-Si [116]. The generally accepted reason for the cube-like shape is the 
effect of anisotropic elasticity. A further observation is the alignment of these 
precipitates, like strings of pearls, along the elastically soft directions. 

As an example, we show in Fig. 1 typical transmission electron microscopy 
(TEM) and small-angle X-ray scattering (SAXS) data for Ni-Al-Mo alloys. This 
ternary alloy system has the property that the lattice misfit depends on the 
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molybdenum content in such a way that negative, positive or zero misfit is 



possible [105 . Fig. la shows the TEM data for an alloy with no misfit. The 
7' precipitates, which appear white, are round in this case, and arranged at 
random. The SAXS pattern (Fig. lb), which represents roughly the square of 
the Fourier transform of the microstructure (as depicted in Figla), is completely 
isotropic. This symmetry shows that there is no preferred direction in the 
configuration of the particles. This situation is exactly what one expects for a 
precipitate microstructure which coarsens to reduce the total amount of interface 
between matrix and precipitates. 

For comparison, Figs, lc and Id show the corresponding data for a Ni-Al-Mo 
alloy with similar fraction of 7' phase but where the molybdenum content was 
adjusted to give a lattice misfit. Now the precipitates are cube-like rather than 
spherical and they are arranged like strings of pearls along the cubic directions 
[100] or equivalent ones. Correspondingly, the SAXS patterns are no longer 
spherically symmetric, although they still have cubic symmetry. Fig. le and Fig. 
If show a similar situation but with a higher volume fraction of 7' precipitates. 
Strongly periodic arrangements of precipitates are also found in alloys like Ti-Ni 

In some other types of alloys, a typical example being Al-Cu jl7] it is also 
common to find plate-like coherent precipitates as a result of elastic misfit 
strains. The example of the so-called 'Guinier-Preston zones' in Al-Cu is par- 
ticularly striking, since these coherent precipitates are practically monoatomic 
layers of copper which create large elastic distortions in the aluminium ma- 
trix. Other examples of alloys with plate-like precipitates are Cu-Be |38| and 
internally oxidized Cu-Fe [f57l p3[. 



A striking effect has been observed in alloys of the type Co-Pt |88 or (CuAu) 



Pt |47fl where precipitates with tetragonal structure are developing in a cubic 
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matrix. Fig. 2a (from p|) shows the tile- like structures that can develop. 

Finally, it has been found in some nickel-base alloys that large cube-like 



precipitates may even split into several smaller ones |59|, 143 1 a process which 



has also been related to elastic interactions [ 1 24 1 (see Fig. 4) 



1.2.2 Effects on coarsening kinetics 

The reported effects on the rate of coarsening are seemingly less consistent. In- 
deed, an increase of the mean precipitate radius proportional to i 1 / 3 , t being 
the annealing time, has been reported in some cases, like Ni-Al-Mo |)7], ||(J or 
Ni-Al-Si |5(|. This growth law is very similar to what is usually found in alloys 
without lattice misfits and in computer simulations. It is due to a ripening pro- 



cess, described in the Lifshitz-Slyozov- Wagner (LSW) theory [198, 199 1, where 
the large precipitates grow at the expense of the small ones, thereby reducing 
the total interface area between matrix and precipitates. 

In some other cases, e.g. Ni-Cu-Si alloys, the coarsening starts out according 
to a i 1 / 3 but slows down in later stages, eventually even coming to a standstill. 
This behaviour is illustrated in Fig. 3 where the growth of the mean precipi- 
tate size R is plotted for Ni-Al-Mo (Fig. 3a) and for Ni-Cu-Si (Fig. 3b) alloys. 
Another example of this type are titanium-rich w-phase precipitates in Ti-Mo 
[ [55l [87f (see Fig. 3c). The usual interpretation of this slowing down is inverse 
coarsening, which means that, because of the elastic interaction, smaller pre- 
cipitates may grow at the expense of larger ones even though this increases the 
interfactial area. 

There is, however, no generally accepted explanation for the fact that slowing 
down of coarsening is observed in some cases but not in all. One possibility is 
that the inverse coarsening may reduce the prefactor A in the expression 

i?«Ai 1/3 (1) 
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instead of affecting the exponent 1/3. This is, indeed, observed in some Ni-Al-Si 
alloys |5(| Another proposal is that there could be an interplay between elastic 
heterogeneity and anisotropy such that inverse coarsening would only become 
predominant in cases where the anisotropy is not too large Q, pj| . 

1.2.3 The role of applied stress; rafting 

When an external stress is applied to an alloy it can cause the shapes and 
arrangement of the precipitates to change. For example, as already mentioned 
and shown in Fig. lc-f, in a nickel-base alloy annealed without external stress 
the precipitates are approximately cubes and are arranged in a cubic array, in 
conformity with the cubic symmetry of the crystal. If a tension is applied along 
the (100) axis the precipitates shorten along this axis and widen along the other 
two axes. At large enough tension the precipitates join up and form plates at 
right angles to the (100) direction. On the other hand a compression along the 
(100) axis causes the precipitates to lengthen along this axis and eventually form 
rods. This behaviour is shown for the case of the Ni-Al-Mo alloy system in Fig. 
lg. The corresponding SAXS pattern (Fig. lh) shows clearly the breaking of 
the cubic symmetry. Pictures very similar to Fig. lh have been published very 
recently for a uniaxially strained commercial nickel-base superalloy studied by 
small-angle neutron scattering 

This phenomenon, known as rafting, was first observed experimentally by 



Webster and Sullivan |186|, Sullivan et al. |183[ and Tien and Copley |176|. A 
well-written short review of the subject is given in the paper by Nabarro et al. 
[ p8| . The mechanism of rafting is thought by many authors (Socrate and Parks 
1993 [||, Pollock and Argon 1994 @, Valles and Arrell 1994 ||, |§, Bufliere 
and Ingot 1995 Q, Svoboda and Lucac |fl| Ohasha et al. 1997 |l§) to involve 
not only the diffusional and elastic effects considered here but also dislocations 
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and plasticity; however such considerations are beyond the scope of this paper. 

2 The static sharp interface model 

In the sharp interface model the two phases, which we call a and /?, are modelled 
as distinct regions of space which, together with the sharp interface between 
them, make up the entire region fl occupied by the specimen. The two regions 
will be called fl a and f?' 3 , and the interface between them will be denoted by T. 
It frequently happens that one of the phases consists of disconnected domains 
called precipitates or inclusions, while the other phase, called the matrix, occu- 
pies the remaining region which looks something like a Swiss cheese. In this case 
we shall normally take Q a to comprise the inclusions and Q.P to be the matrix. 

There are two versions of the sharp-interface model, which we shall call 
the static and dynamic versions. In the static version, time evolution is not 
considered, and the local composition of the alloy need not be considered either; 
the main thing that the model can do is to compare the (free) energies of various 
configurations of inclusions. In the dynamic version, which we shall consider in 
the next chapeter, we do consider time evolution, and since the mechanism 
for time evolution is diffusion the local composition of the alloy must also be 
included in the model. 

Although the static model contains no explicit time evolution and therefore 
no coarsening, it does provide information about time evolution and coarsening 
through the principle that the free energy must decrease with time. For example, 
we can look for the energy-minimizing shape of an isolated inclusion of a given 
volume, and under suitable conditions the actual inclusions having this volume 
should be close to this shape (which is not necessarily spherical, even in an 
isotropic material). The theory of the model also contains important general 
results about the elastic energy, such as the 'Bitter-Crum' theorem (see section 
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|2.1.4 ), and the elasticity theory contained in it underpins the models described 
in later sections. 

2.1 Theory of the energy of inclusions in an elastic medium 

The free energy of the system can be written in the form 



where / denotes the thermodynamic free energy density at zero stress w is 
the elastic free energy density (defined by making f + w the total free energy 
density), F r denotes the free energy of the interface and W ext denotes the 
potential energy of any external mechanism that the specimen may be connected 
to, for example a weight attached to one end of the specimen. In the static 
version of the model, which we are considering in this chapter, the value of / 
can be different in the two phases, but is uniform throughout each phase. The 
variable of integration x is the coordinate of the relevant material point in the 
unstrained reference state of the material, so that a given material point always 
has the same value of x, regardless of where the deformation has taken it to. 

To write down a formula for the elastic energy density w we first define the 
displacement field u; its definition is that the material point x in the unde- 
formed material]^] will be found at the point x + u(x) in the deformed material. 

2 There is some arbitrariness in the choice of an 'undeformed' lattice, since a uniform 
expansion of an undeformed lattice gives another undeformed lattice. Different choices of 
undeformed lattice are related by a dilatation of the coordinate system, which adds the same 
constant multiple of the unit tensor to dui/dxj and to the stress-free strain tensors. This 
transformation leaves the energy expression ^ invariant or (if the externally applied stress, 
modelled as in eqn includes an isotropic component) merely adds to it a constant which 
is independent of the state of the specimen. Physically measurable quantities such as the 
energy differences between different configurations and the differences of the stress-free strains 
between the two phases are not affected by the choice of undeformed lattice. 
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The strain tensor is then defined (we are using linear elasticity theory; for the 
nonlinear formula see eqn (1,3) of J20l| ) by giving its Cartesian components^] 

, . If dui dui \ 
^I^M + t 1 (3) 



k dxj dxi 

We shall also need the stress- free strain tensor Q e?-(x), defined as the value 
which the strain tensor at x would have if the material were uniform and un- 
stressed. The stress-free strain is taken to be a different constant in each of 
the two phases, corresponding to their different chemical compositions. The 
difference between the stress-free strains in the two phases is called the trans- 
formation strain or misfit strain. A common assumption is that, under zero 
stress, one of the two phases is simply a scaled-up version of the other. Let us 
denote the scale factor by 1 + q, so that the stress-free lattice spacing of the 
inclusion material is greater than that of the matrix material by a factor 1 + q, 
corresponding to a relative difference in specific volume of approximately 3q if 
q is small. Then the stress-free strain can be given the value zero in the matrix 
(this corresponds to a particular choice for the reference state of the material) 
and to have the simple form 

4 = ( 4 ) 

in the inclusions. 

In linear elasticity the elastic energy is a quadratic function of the strain 
tensor. If we adopt the convention that the elastic energy is zero at zero stress, 
we can write this quadratic function 

W = \Y1 A yro«( e w( X ) ~ 4( X )X e ™»( X ) " e m«( x )) ( 5 ) 
ijmn 



3 Some authors define e^j for i j without the factor 1/2 in eqn (|3|), but if that is done 
is not a tensor. 

4 This tensor has also been called the 'intrinsic strain' [ L4C | , and the 'spontaneous defor- 



mation' [ 137 1 . 



15 



Here Xijmn, a fourth-rank tensor, is called the (isothermal) elasticity tensor or 
stiffness tensor. This tensor is positive definite, in the sense that w is positive 
unless the tensors and are equal. In general the components Xij mn may 
vary with position (and with chemical composition, temperature, etc), but in 
the sharp-interface model they are taken to be constant in each phase, though 
not necessarily the same constant. 

The Xijmn $ are symmetric under interchange of i and j, of m and n, and 
of the pair with the pair (m,n); and they usually have further symmetries 
reflecting the symmetries of the crystal lattice. The most symmetrical case of 
all is an isotropic solid such as rubber, for which the elastic energy density 
simplifies (in linear elasticity theory) to 

w = ^{E( Ae W 2 + G Ei( Ae )^ - l s v E( Ae )^) 2 ( 6 ) 

k ij k 

where K is the bulk modulus, (Ae)jj means — and G is the shear modulus 
(often denoted by fi) . The non-zero elements of the elasticity tensor in this case 
are Ann = K + |G, A1122 = K — |G, A2323 = G and the ones related to these 
by symmetry, e.g. A2332 = Ai2i2- A more realistic case for metallurgists is the 
cubic crystal, for which in the conventional notation 

w = ^C*iiE( Ae ) 2 4 + c 'i2{(Ae)ii(Ae)22 + (Ae)22(Ae)33 + (Ae)33(Ae)ii} 

i 

+ 2C 44 ((Ae)? 2 + (Ae) 2 23 + (Ae&) (7) 

so that Ann = C n , A1122 = Ci 2 , A2323 = G44, etc. 

The formula (1q), which generalizes a standard formula of elasticity theory 



(eqn (10,1) of [201]) to the case of non-vanishing stress-free strains, is due to 
Khachaturyan [ 188 1 and (independently) Roitburd |185| |. The formula corre- 
sponds to a particular way of apportioning the total free energy density between 
the two terms / and w, in which / is defined to be the free energy at zero stress. 
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The alternative convention, to define / as the free energy at zero strain, has 



also been used [ 101 1 . 

For the interfacial energy term F r in (||) we shall use the simplest reasonable 
assumption, 

F r = <r\r\ (8) 

where a is the surface tension, assumed constant, and |T| = J r d 2 x denotes 
the area of the undeformed interface. More general possibilities have also been 



considered — anisotropic surface energy by Cahn and Hoffman (1974) [ 168 1, 
surface energy depending on the deformation as well as the undeformed area by 



Gurtin and Murdoch (1975) p6| , Cahn and Larche (1982) JT5lj — but these 
refinements have not been important in the theory of coarsening. 

To model the externally applied forces in a simple way, we shall assume 
that there are no body forces, and that the tractions (forces per unit area) 
applied to the surface <9f2 of the specimen are of the form V . £?f*n,-, where n, 
is the unit outward normal vector and ifj , which we assume to be symmetric 
and independent of position, is the externally applied stress tensor. The energy 
W ext of the mechanism providing these tractions can then be equated to its 
energy when the specimen is in the reference state minus the work it does in 
bringing the specimen from its reference state to the state with displacement 
field Ui, so that 



W ext = const. 



/ ]>>^;Sd 2 x (9) 



Choosing the constant to be zero, and using first the divergence theorem (applied 
to the vector field obtained by multiplying Ui by a unit vector in the j direction) 
and then the formula (|^) and the symmetry of tf xt we can write this in the 
alternative form 
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E^ d3x - ( 10 ) 



n « 



In this formula, is the strain in the presence of the externally applied stress, 
which in linear elasticity is an affine function ^| of t\ xt . With neglect of terms 
quadratic in the externally applied stress, therefore, the strain in ( |l0| ) can 
be given the value it has when there is no externally applied stress. 

2.1.1 Elastic equilibrium 

To evaluate the elastic terms in (||) , we need the displacement Ui at each point 
in space. The condition determining itj is that it should minimize the elastic 
energy J fi u>d 3 x + W ext . The Euler-Lagrange equation for this minimization, 



called the equation of elastic equilibrium (eqn (2,6) of [201]), is (since tfj is a 
constant) 

E^=o in, 

3 J 



where t%j is the stress tensor, defined (eqn (3,6) of 201 ) by 



9w , . 

* = (12) 

and Uij means dui/dxj. 

Using (|5|) in ( |l2| ) we find with the help of the symmetry relation Xij mn = 
Xjimn that the stress tensor is related to the strain tensor by 

mn 

(Hooke's law) and is therefore symmetric. In the isotropic case, with w given 
by (||), Hooke's law takes the form 

tn = K5a ^(Ae) fcfc + 2G{(Ae) y - ^(Ae) fcfc } (14) 

k k 



5 That is, a constant plus a linear function 
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To complete the specification of the elastic equilibrium problem, we need 
conditions on the elastic field variables at the interface T. and at the boundary 
of the specimen. On T, the condition of coherency requires u to be continuous: 

[u] = (15) 

where [u] means [u]^, that is u" — u' 3 , the difference between the limiting 
values of u on the two sides of the interface. Moreover, the values taken by 
the field u at the interface (which describe the position of the interface in the 
deformed material) must satisfy an energy minimization condition. Since we are 
assuming that the energy of the interface depends only on its position relative 
to the lattice, the surface energy is not changed by changes of the value of u at 
the interface, and the condition of minimum elastic energy turns out to be 

5>>, = (16) 

i 

where [tij] means tfj — A, and n is a vector normal to the interface; the physical 
interpretation of eqn ( |l6| ) is that the tractions exerted on the interface by the 
material on each side of it must be equal and opposite 0. 

The boundary conditions at the surface of the specimen depend on the physi- 
cal situation there. If the boundary of the specimen is clamped, then the system 
of equations (^Tj) and ( |l3| ) should be solved subject to prescribed values of u. k 
on dfl. In the more common situation where known forces are applied at the 
boundary the boundary condition (obtained by minimizing F without restrict- 
ing Ui at the boundary) is that the surface tractions due to the internal stresses 
are equal to the externally applied surface tractions : 

X ( 17 ) 



If the interfacial energy depends on the deformations as well, then an extra term appears 



in (hfl); see Gurtin and Murdoch (1975) [166], Cahn and Larche (1982) [151 
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where, as before, rij denotes a unit vector normal to the boundary <9f2. 

It is also possible to use periodic boundary conditions, taking Q to be a cube 
(or rectangular prism) and requiring the stress and strain fields to be periodic 
with unit cell fl. In general, the displacement field m will not be periodic, but 
by virtue of (^]) it can be written in the form 

u t = u\ + ^ UijXj (18) 

3 

where v! i is periodic ^| and C/y is a constant tensor, not necessarily symmetric. 
Then, using (|^), the strain can be written in the manifestly periodic form 

(i9, 



p / e 4J (x)d 3 x^(l^ +[/,,) (20) 



2 \dxj dxi 
where e»j is the space-averaged strain tensor 

1 f, — 1 
In 

The value of does not appear in the elastic equilibrium equations but it can 
be determined by direct minimization of F with respect to e^-; the minimum 
occurs when 

[ t (J d 3 x= f t\f d 3 x (21) 
Jn Jn 

2.1.2 The interface condition for inclusions of fixed volume 

The position of the surface T (relative to the undeformed material) is deter- 
mined by the way the atoms diffuse through the specimen. As they diffuse, the 
position of T will change, but the total volume of each of the regions £l a and 
Q.P will remain the same. In the final equilibrium configuration, the position 
of r will minimize the energy subject to given values of these two volumes. In 
addition, if the various inclusions are far enough apart, the volumes of the in- 
dividual inclusions will change only slowly, and it is a reasonable medium-term 



7 Khachaturyan [15C] calls u£ the 'local displacement' 
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approximation to minimize the energy with respect to T keeping the energy of 
each separate inclusion fixed (assuming, of course that such a minimum exists). 
In this section we formulate a necessary condition which T must satisfy at any 
such energy minimum. 

Let r be the surface of an inclusion f2 a in a matrix fi' 3 and consider a small 
displacement of T to a new position V, each point x on T moving to point 
x + n(x)fe(x) on V (see Fig. 5). Here n(x) is a unit vector perpendicular to 
r at x, pointing outwards from the inclusion (i.e. from £l a to Vl^ , and Sx(x), 
a small scalar, is the signed perpendicular distance from x to the nearest point 
on r', positive if T' is in the (3 phase of the T configuration. Denote the region 
between T and T' by SQ a . Since we are not considering diffusion at present, the 
volume of the inclusion is fixed, so that the volume of Sfl a , with an appropriate 
sign convention, is zero. Thus we have (to lowest order in Sx) 

J 5xd 2 x = (22) 

Now consider the change in energy on going from r to V . There is no change in 
the total thermodynamic free energy J f d 3 x, since the volume of the inclusion 
does not change. The change (increase) in the elastic energy between the two 
configurations, as given by eqn (j^), is the sum of three terms: 

1. the change in surface energy due to the change in the area of T. This 
energy increase is equal (in three dimensions) to the surface tension a 
times the increase in undeformed area J t (2kSx) d 2 x, where k is the mean 
curvature, taken positive if the centre of curvature is on the a side of T. 

2. the change in total elastic energy due to the fact that (for positive 8x) the 
region (5f2 Q , which was previously on the (3 side of the interface, is now 
on the a side. To lowest order in Sx this increase is J r [w]5xd 2 'x. where 
[w] = w a — w 13 is the difference between the elastic energy densities on the 
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two sides of the interface at x. 

3. The change in the elastic energy due to changes Su in the displacement 
field outside Sil a . Denoting these displacements by Su and using the 
definition ( |l2| ) of we can write the energy change, to lowest order, as 



Ujj d 3 x 



= / /.{(tjjSuj)^ - tjjjSuj} d 3 x (23) 

where the notation (. . .), j means d(...)/dxj . The second term in the inte- 
grand is zero because of the elastic equilibrium condition (p~T|) . Applying 
the divergence theorem separately to the regions n a and nP, and assum- 
ing that the surface of the specimen is either held fixed, so that Su = 
there, or is too far away to matter, we can simplify the remaining term to 

/ 5>iM«i]d 2 x = f ^n^lSuijd 2 ^ (24) 

by (|l6|). The last equation would also be true with replacing tfj. To 
evaluate [Sui] (which means we write it as a line integral across 

Sfl a 

[Su i ] = [Su i }% = -[Su i f a = J^Y^^-^dz 

Jo ^ ' '•' . 

= — Sx2^[iii,j] n j (25) 

j 

Adding together the three contributions we find for the change in the total 
energy W 

SW = j {2<jk+ [w] - VT 4 [u lfc ]n fc }fed 2 x (26) 



dxj J 



22 



where we have defined Ti = V\ t%rij = V . tf-n,-. the traction at the interface. 

For W to be a minimum against arbitrary changes of shape or position of 
the inclusion which satisfy the constant-volume condition ( p2] ) it is necessary 
for the coefficient of 8x in the integrand of (|2^) to satisfy 

2<jk + [w] - Ti[u i; k]nk = P (27) 

where p is a Lagrange multiplier, the same at all points on T, which can be 
thought of as the excess pressure inside the inclusion. If there is more than one 
inclusion, p can be different for different inclusions. 

If the inclusion is not in equilibrium with respect to changes of shape or 
position, then its surface will tend to move, since the energy can be decreased 
by moving the interface towards the a phase in the places where the left side 
is greater than the right side, and towards the (3 phase in the places where 
the right side is greater. The mechanism for this movement of the interface is 
diffusion, which we discuss in chapter |[ 

2.1.3 The Fourier transform solution 

If the elasticity tensor is independent of position (i.e. , for the sharp interface 
model, if it is the same in both phases) we speak of homogeneous elasticity. In 
that case, the equation of elastic equilibrium, ([n]), can be solved by Fourier 
transforms for the case where f2 is a rectangular box with periodic boundary 
conditions. In this way the elastic energy of an arbitrary arrangement of in- 
clusions in such a box, or an arbitrary periodic arrangement of inclusions in 
an infinite specimen, can be expressed as a sum or integral in Fourier space. 
Moreover, instead of being restricted to two different values, one in f2" and 
the other in , the stress-free strain tensor e^- can vary with position in any 
way. This method of solution has been much exploited by Khachaturyan; see 
his book p0|. 
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Wc define the Fourier components of the displacement vector and the stress- 
free strain tensor by 

u(k) = f e ik x u'(x) d 3 x (28) 

e° mn (k) = [ e ik x e° (x)d 3 x (29) 

where u[ is defined in (|l8|). Using Parseval's theorem and the symmetries such 
as Xijmn — Xjimni the elastic energy formula implied by (|^) and (|l(i| ) can be 
written 



W( eij (-)) = [ {w-Y.nfeij)^ 
Jo. 



= iE .E K*( k )^ - 4(k))A ijro „(-i^(k)fc Ti - e^„(k)*) + 

' ' ky^O ijmn 

+ ^ T,(^-^ X ^(e mn -e° mn )-\n\J2t e l fe lJ (30) 

ijmn ij 

where |f2| means the volume of fi, the star denotes a complex conjugate, the 
k-summation goes over the reciprocal lattice of f2, means the space average 
of e.y(x) (see eqn (pp|)) and e?- the space average of e^. Minimizing the right 
side of ( |30| ) with respect to u(k) and ey, we obtain for the actual elastic energy 



(see [|18§ |181| [150 



w = ^EE e°-(k)^ ijm „(k)e^(k) - \n\J2tT (4 + 5 e lf) (31) 

k^O ijmn ij 

where 

(MO), (32) 

Z^- (k) being the inverse of the matrix 

)u ^ fern ^imnj ^n ; (33) 

and e??* means the solution of J2 mn Ayron e mn = ^mn- Note tnat *ijmn is ho- 
mogeneous of degree zero in k : it depends only on the direction of the vector 
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k, not on its length. Also, \&(k) is positive semi-definite for each value of k, 
in the sense that J2ij mn e ij*ij'mn(k)e mn > for any tensor ey; this property 
follows from the fact that this last sum is the value of the obviously positive 
corresponding term of the sum over k in (|3^) for a particular value of u(k). 

In an alloy where e^ (x) takes just two values, eg in Q a and e?- in ft 13 , the 
definition (p9|) gives 

4(k) = [4Mk) (k^O) (34) 
where [e^-] stands for the transformation strain eg- — and 6*(k) is defined by 

0(k) = f e ikx d 3 x. (35) 

Substitution of ([34j) into the double integral in (^lj) gives Khachaturyan's for- 
mula mm 

W-W ext = i^B(k)|0(k)| 2 (36) 

where we define 

S(k) = ^ [e9.]^. m „[e^J (37) 

ijmn 

Like 4'y TOn (k), the function £> depends only on the direction of the vector k, 
not on its length. 



Using Parseval's theorem, the sum over k in (31) can be expressed as an 
integral over position space : 

W =\ /^E(4( x )-4)^n( x -y)( e m„(y)-4n)d 3 xd 3 y+^ ext (38) 



ij mn 



where Hij mn stands for |0| 1 J2^o *ymn(k)f 
form of ^ijmn, and W ext is equal to the sum involving tfj* in (|30|). Eqn ( j38| ) 
shows how, in the case of homogeneous elasticity, the effect of elastic forces is 



equivalent to a certain two-point interaction [15C, 1 19 1 . It is a singular inter- 
action, however : it generally includes a delta-function at zero range, arising 
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from the fact that the k-space average of ^ij mn (k) is not zero. For anisotropic 
elasticity this interaction also decays very slowly at large distances, because its 
Fourier transform depends only on the direction of k and is therefore discontin- 
uous as k — + 0. In three dimensions, iJy m „(x — y) decays like the cube of the 
distance x — y|, for distances small compared to the size of Q. 

In the special case where the material is isotropic (see eqn (^J)) as well as 
being elastically homogeneous and where the stress-free strain is a pure dilata- 
tion as in eqn (Q), only the delta- function component of Hij mn is present, and 
eqn ( |3^ ) simplifies to 

18KG 



W 



3K + 4G 



f (?(x) - q) 2 d 3 x + W ext (39) 
Jn 



where g(x) is defined to be q inside the inclusions and outside, and q is the 
space average of g(x). The result (|39|), without the term W ext , is due to Cahn 



(1961) [ 1 95 1 . He writes the expression in front of the integral sign more neatly 
as 2E/(1 — v), where E is Young's modulus and v is Poisson's ratio.^j 

The question of what formula replaces (B8h if the elastic stiffness matrix 



is heterogeneous (i.e. different in the two phases) was considered by Onuki 



in 1989 119 for isotropic elasticity and in 1990 [102, 101] for a cubic crystal. 
Treating the inhomogeneity as a small perturbation, he found that to first order 
of perturbation the elastic energy could be expressed as the sum of three parts 

1. a long-range interaction due to the anisotropy 

2. a dipolar interaction which is proportional to the external stress and to 
the amount of inhomogeneity 

3. a non-quadratic 'Eshelby' interaction proportional to the amount of inho- 
mogeneity 



8 Young's modulus and Poisson's ratio are related to the shear and bulk moduli by K 



£7/3(1 - 2u),G = £7/2(1 + 1/). See, for example, Gurtin \174\ 
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Khachaturyan et al. (1995) |^| carried this type of analysis further, considering 
higher orders of perturbation. They found that in the lowest order of pertur- 
bation the inhomogeneity is equivalent to a four-point interaction; in the next 
order, to a six-point interaction, and so on. 



2.1.4 Two theorems 



Theorem 1. (Eshelby 1961) [196]). If the elasticity tensor Xij m n and the external 



stress tensor tf?* are independent of position, W ext depends only on the total 
volume of the inclusions and not on their sizes, shapes, number or arrangement. 

Proof By the conditions of the theorem, the tensor efj* defined just after 
eqn ( |33| ) is independent of position. Hence eqn ( |Io| ) can be written 

ij ^ mn 

= - E e if I (*« + E Vnel) d 3 x by (13) 

ij ^ mn 

= - E e if I / d ' X + E ^mn(|fi°K n + \ (40) 

by (|2l|) . This last expression depends on the inclusions only through their total 
volumcsl, QED. 

Theorem 2. (The 'Bitter-Crum' theorem) [] Under the conditions of theo- 
rem 1, if in addition the elasticity tensor and the stress-free strain tensor are 
isotropic, then for an infinite or periodic specimen the elastic energy W depends 
only on the total volume of the inclusions, not on their sizes, shapes, positions 
etc. 

9 It was first pointed out by Bitter [ ^0^ that in in a linear elastically isotropic material there 



is no elastic interaction between two inclusions of a certain type, and Crum [ 207 ] (also, inde- 



pendently, Robinson [206 ) showed how to generalize this result to a wider class of inclusions. 



See also Eshelby 1957, 1959 [205, 200 1. The name 'Bitter-Crum' theorem was, apparently, 



devised by John Cahn [142 
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Proof : For periodic boundary condition, the result follows from the formula 
(f{9|), since the integral is equal to \il a \(q a — q) 2 + \Q,P\(qP — q) 2 and is therefore 
independent of how the inclusions are arranged. For an infinite system, the 
integral in ( |39| ) is replaced by one over all space, and since q is now equal to 
qP, the integral is equal to |f2 Q |(g" — q 13 ) 2 and is again independent of how the 
inclusions are arranged, QED. 

Both theorems, especially the first, are somewhat counter-intuitive. The first 
theorem tells us that in order to explain rafting using linear elasticity theory we 
must either use a heterogeneous stiffness tensor or else have more than one type 
of inclusion so that the stress-free strain can be different in different inclusions. 
The second tells us, as John Cahn has always seen very clearly |3jJ, that if there 
is to be an elastic interaction affecting the shape or separation of inclusions then 

at least one of the following conditions must be violated: 

1. Matrix and inclusion(s) have the same elastic stiffness tensor 

2. The strain characterizing the misfit is isotropic (i.e. purely dilatational) , 

3. The elasticity tensor is isotropic, 

4. The crystal has no boundaries (i.e. it is infinite, or periodic). 

5. The stress depends linearly on the displacement field 

As we shall see, different violations produce different effects on the coarsening 

behaviour. 

2.2 Energy calculations for isotropic elasticity 

In the rest of this chapter we consider some applications of the static sharp 
interface model. In the present section we shall use the approximation that 
both the stiffness tensor and the stress-free strain are isotropic. Although real 
metals are always anisotropic, such calculations are useful because they can 
be used to investigate the effect of elastic heterogeneity (i.e. of violations of 
condition I of the Bitter-Crum theorem). Effects due to anisotropy will be 
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considered in the following section, E 



2.2.1 Elementary calculations with isotropic elasticity and misfit 

We consider first two types of inclusion whose elastic energies in an infinite 
matrix can be calculated quite easily, the plate (slab) and the sphere. To allow 
for the possibility of finite volume fractions, in which case it is not obvious a 
priori which phase will constitute the inclusions and which the matrix, we shall 
allow the stress-free strain to differ from zero in both phases, being given by 
the isotropic formula (|j) but with different values for q in the two phases. The 
object of the calculations is to find the space-averaged energy per unit volume 
and also the space-averaged strain (which gives information about the response 
to an externally applied stress). 
Plate-shaped inclusions 

Consider a configuration consisting of alternating slabs of the two phases, 
the slabs of phase a having thickness h a , those of phase /3 having thickness 
br . Denote the normal to the faces of the slabs by n. The equations of elastic 
equilibrium have a solution in which the strain in the a phase has the form 

e% = Q a $ij + P a {5 lJ - 3n inj ) (41) 

where P a and Q a are constants, and similarly for the p phase. The total energy 
per unit volume of material can be calculated from eqn (||); for an isotropic 
stress-free strain, as given by eqn(|J) generalized to allow for non-zero values of 
q in each of the two phases, it is 

w = {^K a {Q a ~q a ) 2 +6G a {P a ) 2 }0+{^KP(Ql 3 -qP) 2 +6G fJ (pP) 2 }(l-0) (42) 

where <j) — h a /(h a + h 13 ) is the volume fraction of phase a. 

The expression (|42]) is to be minimized subject to the coherence constraint 
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(y_5|), which because of (|3j) and (|41j) is now 

P a + Q a = P 13 + Q 13 (43) 

(This minimization automatically takes care of the stress continuity condition 
(p^l)). The minimum value, which gives the actual mean energy density of this 
system of inclusions, is 

ffl = ?fc« K , [# (44) 

where is defined by 

1 1-0 3/1-0 ) 
X, K« 4( G« G^J' 1 j 

In the limit where the volume fraction of the precipitates goes to zero, (this 

includes the case of a single plate-like inclusion in an infinite matrix) the energy 

per unit volume of precipitate is 

. , 18G a K a [q} 2 
lim w 4> = — (46) 

If there is a small externally applied stress, the resulting energy can be 
found from the average strain at zero applied stress, using (|l0|). The anisotropic 
(deviatoric) part of this average strain is 

*a - \ E ^ = -^ pa + o- - ^ = \(&-w) ^ 

k ^ ^ 

We shall not need the isotropic part. 

It is not hard to generalize these results to a cubic crystal, if n is in a 
convenient direction. If n is along one of the crystal axes, it is only necessary to 
replace K by (Cn + 2Ci2)/3 and G by (Cn — G\%)l2 in the various formulas. 
If n is along the [111] direction, G should instead be replaced by C44. 
Spherical inclusion 

Consider a spherical inclusion of a phase, of radius R\, embedded in a larger 
sphere of j3 phase having radius Putting the centres of both spheres at the 
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origin, we look for a solution in the form 

Q«X (|x|<i?!) 

(48) 

Q^x + PV^ (J? x | <x| <R 2 ) 

where Q a ,Q p and P are scalar constants. Since u is radial in both solutions, 
the coherency condition ([Tq) can be satisfied by making 



Q a Rx = QPRx - P/Rf (49) 
Using (|J) we find that the strain is 

Q a 3ij (|x| < Ri) 

(50) 

+ PfaiXj - |x| 2 <5y)/|x| 5 (i?i < |x| < R 2 ) 
Hence, by (^]), the elastic energy is 

4W?i f 9^°= ^ , . 9KP. 

3 

where ^ = (R1/R2) 3 is the volume fraction of phase a, and we have used ( [f9|) 
to eliminate P in the second line . 

The expression ( |5l| ) is to be minimized with respect to Q a and Q* 3 . The 
minimum value, which gives the actual space-averaged elastic energy per unit 
volume, can be written 

^««.H' (52) 

where K m is defined by 

1 * 3 (53) 



{^f(Q a - i a ) 2 ^ + —7j-(Q p - /) 2 (i -<!>) + 6G^(Q' 3 - Q Q ) 2 <Mi - <^x}i) 



K. K a KP AG?' 
In the limit of zero 0, which includes the case of a single spherical inclusion in 
an infinite matrix, the energy per unit volume of inclusion is 
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For this geometry, the space average of the strain is purely isotropic, and so 
the linear part of the energy of a spherical inclusion due to an externally applied 
field is zero. 
Conclusions 

^From these elementary calculations we can draw some tentative conclusions 
about the effects of elastic inhomogeneity on the energetically preferred config- 
urations of the precipitates, provided that that elastic anisotropy does not mask 
these effects. 

First of all, a comparison of eqns ([34]) and ( fSf ) shows that if G a > G 13 , 
that is if the inclusion is more rigid than the matrix, then the sphere has lower 
energy than the plate; on the other hand if the inclusion is less rigid than the 
matrix, then the plate configuration has the lower energy. That is to say, a 
spherical inclusion is likely to be stable only if it is stronger (more rigid) than 



the surrounding matrix. If it is justified to use eqn (52) to estimate the energy 
of a system of spherical inclusions even when <f> is not vanishingly small then 
one can go even further : this formula, when compared with ( [44] ) , suggests that 
a system of spherical inclusions made from the more rigid material will always 
have a lower energy than the alternating plate configuration with the same <j>, 
so that in addition to destabilizing soft spherical inclusions elastic heterogeneity 
tends to destabilize all plate- like inclusions. 

Secondly, eqn (^) shows that for an alloy where the more rigid material also 
has a larger stress-free volume (the correlation between size and strength we all 
remember from the school playground) the space-averaged strain is a contraction 
in the direction perpendicular to the plate and an expansion in the directions 
parallel to it. Consequently in such a material, a plate-shaped inclusion will tend 
to align perpendicular to an externally applied uniaxial compressive stress, but 
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parallel to an externally applied tensile stress. ^| In symbols 

taxial [<?] [G] > =4* inclusions parallel to axis of applied stress 

tlxlai [l] [G] < =>■ inclusions perpendicular to . . . (55) 

where t^\ al = (2t|3* — tff* — t^)/^ is the axial component of the anisotropic 
part of the applied stress (taking this axis to be along the 3 axis), which is 
positive (negative) for a tensile (compressive) stress. 

These conclusions can be understood physically in the following way: in the 
case of a spherical inclusion, the matrix is sheared but the inclusion is not, so 
the sphere is a good (low-energy) configuration when the matrix is very flexible. 
On the other hand, for a plate-shaped inclusion, only the inclusion is sheared, 
and so when the inclusion is very flexible the plate is a better shape than a 
sphere. Moreover if, in addition to being softer than the matrix, the plate has a 
lower stress-free density, then the matrix will prevent it from contracting in the 
direction parallel to the plate; it can only contract in the direction perpendicular 
to the plate, and hence the anisotropic part of the strain is a contraction in this 
direction, and plate-like inclusions will be energetically favoured when the stress 
is compressive. 

2.2.2 The energy of a single inclusion : isotropic elasticity 

For most shapes of inclusion other than sphere, cylinder and plate the calculation 



of elastic fields is complicated, but Eshelby [202, 20C, 196] showed how to do it 
for an ellipsoidal inclusion. He found that if the inclusion is more (less) rigid 
than the matrix, an ellipsoid has more (less) elastic energy than a sphere of the 
same volume. Thus, a spherical inclusion which is more rigid than the matrix 

10 Thc formulation (|S5|) is given in According to Thompson and Voorhees this 

effect was predicted by Cahn soon after his 1961 work on spinodal decomposition : see ref. 
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is stable against deformation to an ellipsoid, but (disregarding the stabilizing 
effect of the surface energy) one that is is softer than the matrix is unstable. 
This result is consistent with the comparison in the previous section between 
a spherical inclusion and a plate, which can be regarded as an infinitely flat 
ellipsoid. 



Very few authors have used the condition (27) for the equilibrium shape of an 



inclusion, but Roitburd [137 has shown that in the case where the transforma- 
tion strain is not a pure dilatation (thus violating condition 2 of the Bitter-Crum 
theorem as well as condition 1), and where the surface energy is neglected, an 
ellipsoidal inclusion satisfies this condition over a certain range of values of the 
parameters. 

If we take surface energy into account as well as elastic energy, the energy 
per unit volume of inclusion depends on the size of the inclusion as well as its 



shape. Johnson and Cahn (1984) [ 146 1 showed that for small inclusions the 
surface term F r in the free energy formula (^) predominates and the shape of 
the lowest-energy inclusion has the same symmetry as that which minimizes F r 
for a fixed volume of the inclusion. In the case of an isotropic surface tension 
this would be a sphere, but for anisotropic surface tension it could also be a 
polyhedron facetted parallel to the crystal planes where the surface tension is 
minimum. As the size is increased the elastic energy becomes more important 
and there is a transition, which may be either continuous or discontinuous, to a 
new equilibrium shape of lower symmetry. If an external stress is present, this 



too can induce a shape transition (Berkenpas et al. (1986) [132], Johnson et al. 
(1988) (ml). 



With a view to constructing a theory of rafting, Pineau (1976) [165 1 calcu- 
lated the energies of single elliptical inclusions in an externally applied stress 
field, using methods devised by Eshelby |196|| . A brief comparison of his predic- 
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tions with experiment is given by Socrate and Parks (1993) ]6^ |. 



2.2.3 The energy of two or more inclusions : isotropic elasticity 

In the absence of elastic forces, the mechanism of coarsening arises from the 
fact that a pair of spherical inclusions with given total volume have less total 
surface energy the more unequal their sizes. To see how this mechanism is 
influenced by elastic forces, we would like to know how the elastic energy of a 
pair of neighbouring inclusions depends on their size and separation. Can they 
reduce their energy by transferring material from one inclusion to another, or 
by changing their separation? 

For isotropic elasticity, calculations of this type were first done by Eshelby 



in 1966(187 . For a pair of spherical inclusions with radii R%, i?2 and centres a 
distance D apart, he found the energy of interaction to be approximately 

w - = 1 {w^h + j^h + ° ([G1)2 } <56) 

where v = (3K — G)/(6K + 2G) is Poisson's ratio and [q], [G] have the same 
meaning as in eqn (f47|). If the inclusions are more rigid than the matrix, then 
formula (|56| ) says that the interaction energy is positive and (for fixed values 
of the separation D and the total volume (47r/3)(i?f + R2)) is least when the 
two inclusions are of equal size. This suggests that the elastic energy will tend 
to make the inclusions of equal size. This stabilizing effect, just the opposite of 
the tendency of surface energy to enlarge the large inclusions at the expense of 
the small ones, has been given the name 'inverse coarsening' 



In his papers drawing attention to the above stabilizing effect Johnson [149, 



145 1 used a generalized version of the formula (|56| ) which included also surface 



energy and the effect of an external stress. He concluded [|149| that, in the 
absence of external stress, the inverse coarsening effect would occur whenever 
the size of the inclusions exceeded a certain threshold related to the critical size 
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we have called Rq. He also found [145] that a uniaxial external tensile stress 



could favour or inhibit inverse coarsening depending on the angle between the 
separation of the inclusions and the axis of the external stress. 

The elastic interaction of more complicated systems of inclusions of given 



shape has also been studied; for example Johnson and Voorhees (1987) [130 
consider the elastic interaction energy of a system of several cuboids (rectan- 
gular prisms), finding a tendency of such precipitates to align, even though the 
medium is isotropic. Abinandanan and Johnson (1993) [^3], Q calculated the 
elastic interaction energy of two spherical inclusions with a tetragonal misfit 
strain. The energy is positive when the line joining the particles and the axis 
of the misfit are parallel, negative when they are perpendicular. 

The value of all such conclusions based on the assumption of spherical pre- 
cipitates is, however, called into question by the work of Onuki and Nishimori 
(1991) J90|. They showed that the interaction energy is extremely sensitive to 
the shapes of the inclusions. Eqn ( p^ ) assumes that the inclusions are spheres, 
but a small deformation of the spheres reduces the interaction to zero. Since 
there is no reason why the actual precipitates should be exact spheres, formula 
( |56| ) cannot be relied on. 

Two other calculations using pre-assigned shapes of inclusion are the energy 
calculations for a three-dimensional array of inclusions in the shape of square or 



rectangular plates by Perovic et al. [ 159 1 and the elastic stress calculations for 



two-dimensional arrays of square plates by Glatzel and Feller-Kricpmayer [ 109 1 . 



2.3 Energy calculations for anisotropic elasticity 

In this section we use energy methods to investigate the effect of elastic anisotropy 
(usually with cubic symmetry) on coarsening — that is to say, the effect of vio- 
lating condition 3 rather than condition 1 of the Bitter-Crum theorem. Except 
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where otherwise stated, the material is assumed to be elastically homogeneous. 



2.3.1 The energy of a single inclusion : anisotropic elasticity 



Khachaturyan (1966) JL88[ [150| (see also Roitburd (1967) |85j) showed how 
to calculate the elastic energy of an inclusion of arbitrary shape, using the 
Fourier transform formula (pi]). This energy is minimized by making the in- 
clusion a plate (or an array of parallel plates) with faces (the 'habit planes') 
perpendicular to a vector n that minimizes the function -B(n) defined in (|37|). 
This configuration obviously satisfies the interfacial equilibrium condition ( p7| ) 
too. The directions of the minimizing vectors are called elastically soft direc- 
tions. For a cubic crystal such as tungsten which has positive anisotropy, i.e. 
C\x — C\2 — 2C44 > 0, the elastically soft (minimizing) directions are [111], [111] 
etc.; for negative anisotropy, as in Al, Cu, Fe, Ni, they are [100], [010] and [001]. 
The elastic energy of some other shapes of inclusion can also be found ana- 



lytically, for example spheres, ellipsoids [181, 182], and cubes |162j, |150|. If the 
strain in the inclusion is uniform (i.e. for all these shapes except the cube), the 
energy can also be calculated for a material which is elastically heterogeneous 



as well as anisotropic [164] (or section 8.4 of |150j ). 

When surface energy is taken into account, the minimum-energy shape of 
the inclusion is no longer a plate, but it has still been computed analytically in 



some cases. The two-dimensional calculations of Wen et al. (1981) [154, 150 1 
illustrate how the minimum-energy shape, round when the inclusion is very 
small, can pass through a succession of more and more plate-like shapes as the 
inclusion is made larger. For a crystal with square symmetry, as the size of the 
inclusion is increased there is a transition from minimum-energy shapes with 
fourfold (square) symmetry to ones with twofold symmetry [|Bl[ p9| . (Symmetry- 
breaking transitions of a similar kind occur for isotropic elasticity : see section 
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In three dimensions, the minimum-energy shape is more difficult to compute; 
a common procedure has been simply to compare the energies of various easily 
calculable shapes such as spheres, ellipsoids, cubes, cuboids, tetrahedra and 
octahedra and to assume that the actual shape at a given volume of inclusion 
will be similar to the easily calculated shape having the least energy. Such 
calculations [124, 112, 109, [77| show that the precipitate may pass thorough a 
variety of different minimum-energy shapes as its size is increased, and may even 



split into two or more smaller pieces. As we have mentioned in section 1.2.1 



splitting of this kind has been observed experimentally with 7' precipitates in 



nickel-base and other superalloys [153, 13S, 112, 59j] : see Fig. 4. A review of all 
such energy calculations up to 1992 is given by Johnson and Voorhees [[74], (73| . 

The effect of an externally applied stress on the energy-minimizing shape 
of an ellipsoidal inclusion with elastic forces only, a key factor in the theory of 
rafting, was studied by Chang and Allen (1991) |84|, following earlier theoretical 



work by Tien and Copley (1971) [176|, Pineau (1976) ]165[ (mentioned in section 
|2.2.2j above), Miyazaki et al. 1979 ]l5S| ], Johnson (1987) JT2| (who showed that 
the true equilibrium shape could be approximated by an ellpsoid of revolution), 



and Johnson et al. (1988)) [122]. Using realistic elastic constants, Chang and 
Allen find good agreement with observed shapes and orientations. Their article 
also contains a good critical review of the earlier theoretical work on rafting. 
Nabarro et al. (1996) [^8| consider a cubic crystal with inclusions having the 
shape of tetragonal prisms instead of ellipsoids, and they find that eqn ( p5| ) 
holds, with the shear modulus G replaced by its analogue for a cubic crystal, 
which is (see the paragraph after (piTf)) (Cn — C22V2 if [100] is a soft direction. 

The restriction of working within a pre-selected set of shapes (spheres, 
cuboids, etc.) is avoided in the work of McCormack et al. 1992 |7(| who 
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do an unusual type of two-dimensional finite-element computation in which, as 
the size of the inclusion is increased, its shape follows a path of steepest descent 
with respect to energy, not even requiring the topology to stay the same. As the 
size of the inclusion is increased, the round shape first turns into a square, and 
later on a nucleus appears inside the square (or on the surface if the calculation 
forbids changes of topology) and the square splits into two or four. 

2.3.2 The energy of several inclusions : anisotropic elasticity 

As in the case of isotropic elasticity, we study pairs of inclusions in order to 
understand when inverse coarsening can take place, and if it does what relative 
position of the inclusions gives the least energy. In addition, there is the possi- 
bility which we have already mentioned (see Fig. 4) that a single inclusion may 
be able to reduce its energy by splitting into two (or more) parts. 

Several authors have studied these questions by investigating the energy 



of a pair of spherical inclusions. Johnson and Lee (1979) [157] show that two 
spherical inclusions in an anisotropic matrix can reduce their energy if they come 
close together and align themselves along the 'soft' crystal directions. Miyazaki, 



Doi and others [135, |7^, ^ find in addition that inverse coarsening is possible 



if the particles are close together and the misfit strain and particle size are large 
enough. Shneck et al. (1992) give details of the elastic fields for one and 
two spherical inclusions in a cubic crystal with negative anisotropy. 

Inclusions shaped like cubes and cuboids have also been studied, because of 
their relevance to the observed splitting of such inclusions in 7' precipitates. Doi 



et al. (1984) [ 143 1 showed that a single cuboid may have more elastic energy per 



unit volume of precipitate than an array of eight cuboids, or a pair of parallel 



plates. Khachaturyan et al. (1988) [124] compare the energy of a cube-shaped 



inclusion, in a cubic crystal having negative anisotropy, with that of two half- 
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cubes and with an octet of smaller cubes. Depending on the distances involved, 
either the doublet or the octet can have lower energy than the single cube. 



2.3.3 Positional correlations 

We have already noted that a pair of inclusions of given shape often has the 
least energy if the inclusions are close together and are aligned along an elasti- 
cally soft direction. The same tendency was found in arrays consisting of a large 



number of inclusions, by Khachaturyan 1969 [18C, 15C|, and also by Khachatu 



ryan and Airapetyan 1974 [170|. They found, for example, that in a cubic alloy 
with negative anisotropy the minimum-energy array for a system of spherical 
inclusions is a simple cubic superlattice. 



3 The dynamic sharp-interface model 
3.1 Modelling diffusion 

The dynamic version of the sharp-interface model differs from the static ver- 
sion in taking account of the mechanism of time evolution, which is diffusion. 
Since diffusion arises from concentration gradients, we drop the assumption that 
was made in the static version of the sharp-interface model that the intrinsic 
properties of the alloy are uniform in each phase. The composition of the alloy 
will now depend on position; but the elastic stiffness matrix and the stress-free 
strain tensor will still be assumed uniform within each phase, though in general 
different in different phases. 

For simplicity, let us assume that the alloy consists of just two types of atom, 
A and B, so that the local composition of the alloy can be described by a single 
field, the local concentration of A atoms, which we denote by c. We shall assume 
that the concentration of A atoms in the (3 phase is small so that c is small in 
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the (3 phase. In the a phase, we assume that c is close to its zero-temperature 
equilibrium value in that phase which we denote by Cq. If the a phase is pure 
A, then eg 1 = 1, but it is also possible for Cg to be less than 1; for example in 
an Ni-Al alloy the precipitates have a M3-AI structure so that for this alloy the 
equilibrium concetration of Al atoms in the precipitates is Cq = 1/4 . 

3.1.1 The diffusion equation and its boundary conditions 

Within each phase, the concentration varies with time by diffusion, which we 
can model by the standard diffusion equation 



For simplicity, we shall take the diffusivity D to be uniform in each phase, 
though it may take different values in different phases. 

The motion of the interface is controlled by the net rate of arrival of A 
atoms, in accordance with the standard mass conservation formula associated 



where n is a unit vector normal to T and v n denotes the velocity of T (with 
respect to the undistorted lattice) in the direction of n. 

To use these equations we need two conditions to determine the value of c 
on the two sides of T. Since these conditions are also satisfied at equilibrium, 
we can obtain approximations to them by formulating the two conditions that 
must be satisfied on T at equilibrium. One of these conditions refers to the 
possibility of mass transfer across T, the other to the possibility of motion of T 
relative to the underlying lattice. 

The condition for equilibrium under mass transfer across T is that the chem- 
ical potentials of the two kinds of atoms on the two sides should be equal. De- 
noting the chemical potential of A relative to B (i.e. the excess of the chemical 




(57) 



with (||), 



[Dn ■ Vc] 



(58) 
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potential of A over that of B) by fi, we can write this condition 

fj, a = yP on T (59) 

where \i a (or denotes the limiting value of p, at any point on T as that point 
is approached from within the a (or /?) phase. 

If the system is not in equilibrium, extra terms may appear in ( |59|) depend- 
ing, for example, on the velocity of the interface (corresponding to the kinetic 
undercooling term in the theory of freezing and melting). We shall assume 
here, however, that the deviation from equilibrium is small enough to justify 
neglecting such terms. 

In defining the chemical potentials which appear in (|5^), we must take ac- 
count of the clastic contribution to the free energy, defining it as the derivative of 
the total free energy with respect to particle number at constant strain. How- 
ever, since we are taking the stiffness matrix and the stress-free strain to be 
independent of concentration in each phase, it follows from the formula (j^) for 
the strain energy that the elastic contribution to the free energy does not affect 
the chemical potential. The formula giving the chemical potential (of A relative 
to B) is therefore just the same as it would be in the absence of elastic energy: 

(If temperature variations were of any importance, we would of course write 
(df/dejx, but in a metal it is reasonable to treat all processes as isothermal.) 

To obtain our second condition, we consider the condition for equilibrium 
under displacements of the interface. This can be obtained by requiring the 
total free energy to be a minimum with respect to variations of T subject to the 
constraint that the number of A particles, denoted here by N a , is fixed: 

N a = I cd 3 x = const. (61) 
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The argument closely follows the corresponding one in section 2.1.2j , with the 



constant-volume constraint used there (see eqn (p2|)) replaced by (pi]), and the 
resulting necessary condition for a minimum, analogous to (p7fi, is 



2an + [/ + w] - Tjlu^klrik = n[c] (62) 

ik 

where [c] = c Q — c' 3 is the discontinuity in c across the interface. We have 
denoted the Lagrange multiplier by /i in accordance with the general princple 
that the Lagrange multiplier is equal to the derivative of the minimum value 
with respect to the value of the constraint, which in this case is dF/dN a , i.e. 
the chemical potential \i. At equilibrium, /i, fi a and \jt are all the same, and 
as in the case of eqn (^9|) we shall make the approximation of neglecting any 
non-equilibrium terms in (|62|), proportional to, say, the velocity of F, even when 
the system is not in equilibrium. 



The condition («32|) (without the curvature term) is due to Robin [171| and, 



in greater generality, to Larche and Cahn (1978) 1 16 1 ] . The curvature term, 
which incorporates the Gibbs-Thomson effect, was included by Cahn and Larche 

on 



(1982) [151 1, following earlier work by Gurtin and Murdoch (1975) [166 



surface stresses in solids. 



3.1.2 The generalized Gibbs-Thomson condition 



Eqn (|62|) can be brought to a more convenient form, generalizing the Gibbs- 
Thomson formula. For each phase we introduce (following Larche and Cahn 



1973 [172 1) a thermodynamic potential 7r, the grand canonical pressure, referring 
to the properties of this phase at zero stress. Considered as a function of c, it 
is defined by 

tt(c) = c M (c) - /(c) (63) 
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so that 

c?7r(c) dfi 



The interface condition (162) can now be written 



[tt] = 2a k + [w] - 7i[tt i)fc ]n fc (65) 

ik 

where [tt] means ir a — ir° , that is 7r(c Q ) — tv(ct). 

To obtain convenient approximate formulas for fi and n, consider first their 
values at the type of equilibrium normally considered in thermodynamics, for 
which the interfacial curvature is small and the lattices in the two phases are 
incoherent, so that surface and elastic energies play no part in the equilibrium 
conditions. At such an equilibrium eqn (^9|) shows that there will be a common 
value of [i a and yP^ denote this common value by fi eq . Moreover, eqn ( |65|) 
reduces in this case to [tt] =0; that is to say, ir a and tt" have a common value 
at such an equilibrium. Denote this common value by ir eq , and the values of c 
in the two phases at such an equilibrium by c° q and c% q . 

By Taylor's expansion about this point of incoherent equilibrium (and, in 
the second line, the use of (pi))), we have the following approximations which 
are useful in the a phase: 

»(c a ) = y eq +^(c« q )(c a -c« q )+0(c a -c« q ) 2 (66) 
TT(c a ) = 1Teq +c^(c« q )(c a -c« q ) + 0(c a -c« q ) 2 (67) 

where the prime denotes a derivative with respect to c. There is a similar 
approximation for the (3 phase. It follows from ([36]) and(^) that 

n(c a )n eq + c^(c a ) ~ M (c« )) + 0( M (c«) - fx eq ) 2 (68) 

^From this we can subtract the corresponding formula for tt(c^), and using also 
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the condition (59) in the form /i(c Q ) = /i(c^), we obtain 

[w] = [c eq ] (p - /j, eq ) + 0(fi - i± eq ) 2 (69) 

where /x denotes the common value of fi a and ^ . Combining ( |69| ) and (|66|), 
and using (pq), we get the approximate formula 



A 1 - Meg 



(70) 



e9 M '(c«) [c e9 ] M '(c«) 
The formula for c' 3 — c£j is analogous. 

For the model alloy we are considering here, the chemical potential can 
be expanded in virial-type expansions about the two points Cq and Cq. If we 
take Cq to be zero for simplicity, and in a case such as Ni 3 Al neglect the small 
concentration of A atoms on the sites preferred by B atoms, these expansions 
have the form 

ex P = -R-{ 1 + °( c )} 



kT 



exp^M = -^{l + O(c "-c)} (71) 

where k is Boltzmann's constant and T is the temperature. The reason why 
/i(c) appears with a negative sign in the second equation is that in the a phase 
the small parameter of the virial expansion is Cq — c, the concentration of B 
atoms, and the relative chemical potential of these atoms is — n{c). Using ( fFT| ) 
in d7Q), we obtain the approximate formulas 



kT[c eg ] ' eq kT[c eg ] 



(72) 



in which [tt] stands for the expression on the right side of (p5|). An example 
showing how the differential equation and boundary conditions formulated in 
this section are used is given in the following section 

The formula ([72h, reducing in the absence of elastic effects to the Gibbs- 



Thomson formula, is due to Johnson and Alexander [134] (see also [127, 



114 



45 



p3| . p9|), who used it (in a more accurate version) to calculate c at the surface 
of a spherical inclusion with isotropic elasticity. 

3.2 Sharp-interface evolution calculations 

3.2.1 Growth or shrinkage of an isolated spherical inclusion 



In this section we calculate (following Laraia et al [125, 113|) the rate of growth 
of a spherical inclusion of radius R in an infinite elastically isotropic matrix, 
with the boundary condition that c approaches a prescribed limit c°° far from 
the inclusion. At the surface of the inclusion its value is given by ([72]) • The 
calculations are simple because the geometry makes [w], and therefore c a and 
c°, independent of position on T, so that they can be treated as constants. 

We make the approximation (which appears to be a good one if RdR/dt <C 
D) of neglecting the time derivative in the diffusion equation (^7|), so that it 
reduces to Laplace's equation. The appropriate solution, with origin at the 
centre of the sphere, is 

c a (r< R) 

(73) 

c°° + {cP - c°°)R/r (r > R) 



where r means |x|. Using (58) we find the rate of growth of the inclusion to be 



^ = Vn = - D [dc/cr] = D ( c /3 _ c oo } 1 (74) 
at H 

If there are no elastic forces, then ( |65|) reduces to [tt] = 2a/ R, and so @ with 
c@ taken from (jT^) gives 

dR_ 2aDcj fl i\ 

dt YT[c eq ]R \R* RJ { ' 

where R*, defined by 2a/R* = kT[c eg ](c 00 - c^ q )/c^ q is the so-called critical 
radius, the radius of an inclusion at whose surface the value of c" would be c°°. 
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If elastic forces are present, then eqn KB) can be evaluated using the R2 — ► 00 



limit of the solution for a spherical inclusion outlined in section 2.2.1. The result 
is 

[tt] = 2a /R + l8K a G {q] 2 /{3K a + 4G" 3 ). (76) 



Eqn (75) still holds, with R* now defined by 

2a/R* = kT[ Ce9 ](c°° - c£,)/c£, - l%K a G%f / (3K a + 4G^) (77) 

Eqn ( |75| ) is a central component of the LSW theory of coarsening at very 
low volume fractions of precipitate, which we have already mentioned and shall 



discuss further in section 3.2.4. In that theory each precipitate has its own value 



of R, which changes with time according to eqn (75). The common value of R* 
is determined by requiring the total volume of the precipitates to be conserved. 
LSW argue that both R* and the average radius of the precipitates will grow 
in proportion to i 1 / 3 . Since eqn ( |75| ) is unaffected by the presence of the elas- 
tic term in (|7^) it follows that the elastic term has no effect whatever on the 
coarsening behaviour as predicted by the LSW theory. The elastic effects will 
be revealed only by going beyond the main assumptions made in the above cal- 
culation, which are that the elasticity is isotropic, the precipitates are spherical, 
and their volume fraction is very small. 

3.2.2 A single non-spherical inclusion 



We have already seen in Section [2.2.2 that a spherical inclusion, with isotropic 



elasticity, can be thermodynamically unstable if it is softer than the matrix. If 
the inclusion is growing, a further instability becomes possible, corresponding 
to the one discovered by Mullins and Sekerka |l92| which in the non-elastic 
case appears as soon as a certain critical radius (which depends on the rate of 
growth) is exceeded. This type of instability has been investigated using linear 
theory by several authors 1 103 : 115, 9lJ p5[ . 



47 



The first application of the sharp-interface model to nonlinear time evo- 
lution of the shape of a precipitate was a two-dimensional calculation done 
by Voorhees, McFadden and Johnson (1992) J8C|| , using a boundary integral 
method to solve the diffusion equation in the matrix (they neglected diffusion 
inside the precipitates). They calculated the shape evolution of a single pre- 
cipitate of arbitrary initial shape for an anisotropic system having the elastic 
constants of Ni-Al, choosing the concentration at infinity so that the volume of 
the precipitate remained constant and it approached an equilibrium shape. In 



contrast to the energy calculations described in Section 2.3.1, they did not find 
any symmetry-breaking transitions : the equilibrium shape always had fourfold 
symmetry even if the initial shape did not, and their precipitates never split 
into pieces; but they attribute this to their not having considered such large 
precipitates (i.e. such large values of R/Ro) as the other workers. 

Jou et al. (1997) jl3| did similar calculations, but they used isotropic elastic- 
ity and chose the concentration at infinity so that the precipitate would grow at 
a prescribed rate. The mass influx has a big effect on the shape, favouring the 
formation of dendrites both for relatively hard and relatively soft precipitates, 
whereas the equilibrium shapes are 'squarish'. If there is an external stress field, 
the dendrites perpendicular and parallel to the axis of this field grow at different 
rates. 



3.2.3 Two or three precipitates 

When studying more than one precipitate, or even a single precipitate in three 
dimensions, it is a complicated task to follow the detailed evolution of the bound- 
ary surfaces. To avoid the difficulty, some authors constrain the precipitates to 
be spherical, determining the rate of change of the size and position of each 
sphere from appropriate integrals of eqn (|5^) over the surface of that sphere. 
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Voorhees and Johnson (1988) 1 128 , also Johnson et at (1989) [111], studied 
two spherical precipitates in an isotropic matrix, solving the diffusion equation 
in the matrix using bispherical coordinates. For precipitates that are softer than 
the matrix, they find inverse coarsening if the ratio of the particle radii is less 
than about 2. Diffusion also make the centres of the precipitates move; this 
also happens in the non-elastic case, but here the movement is faster. Johnson 
et al. (1990) |j^] extend this type of calculation to anisotropic elasticity. They 
find inverse coarsening when the preciptate separation is along a soft direction 
([100] in Ni, [110] and [111] in Mo), but coarsening rates are enhanced for some 
other directions of the separation. 

Getting away from the assumption of spherical precipitates, Su and Voorhees 
(1996a, b) |2^, [?0) did 2-D calculations for a pair (or triplet) of inclusions in a 
matrix with cubic symmetry, allowing for shape changes (but not for diffusion 
inside the precipitates). Large precipitates adopted approximately rectangular 
shapes, rather than the round ones assumed in fixed-shape calculations. Inverse 
coarsening was never observed for two-precipitate systems (notwithstanding the 
results of earlier calculations at fixed shape); in systems of three precipitates, 
local inverse coarsening was observed but it was never enough to stabilize the 
system as a whole against further coarsening. With regard to the movement 
of the precipitates, the main effect of the elastic forces was to align pairs of 
precipitates so that the line joining their centres was along one of the elastically 
soft directions. These forces also caused the precipitates to move close together 
(partly by motion of their centres, partly by shape changes), but not to touch. 
In discussing these results the authors make use of the concept of configurational 
forces due to Eshelby ]205[ and Gurtin [0 
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3.2.4 Many precipitates : the LSW theory 

The first theory of coarsening, which took no account of elastic interactions, was 
the LSW (Lifshitz and Slyozov |l98| , Wagner |l99|] ) theory, which is based on eqn 
( [75| ) for the rate of growth of an isolated spherical precipitate. Precipitates that 
are larger than R* in radius grow by diffusion at the expense of those that are 
smaller. The theory, which assumes that the volume fraction of precipitates is 
vanishingly small, predicts that the average radius of the precipitates present at 
time t will grow in proportion to t 1 ^ 3 . For finite volume fractions the LSW theory 
becomes inaccurate, and it can be improved upon by following numerically the 
individual sizes of a system of precipitates which are assumed to be spherical. 



For reviews of this theory see Voorhees (1985, 1992) [141, 



Enomoto and Kawasaki (1988,1989) |L2J, [|K)7| extended this type of nu- 
merical simulation by including the effect of elastic forces. They assumed that 
the precipitates were spherical, so that the total free energy could be expressed 
in terms of the radii of all the inclusions, using Eshelby's formula (^f|) for the 
elastic interaction. (If this interaction is not included, then the clastic energy 



does not significantly affect the rate of coarsening : see [112]). The diffusion 
potential at the surface of an inclusion was replaced by an average value over 
that surface, calculated in terms of the partial derivative of the total free energy 
with respect to the radius of the inclusion. For precipitates that were harder 
(i.e. more rigid) than the matrix, they found that the elastic interaction assisted 
the coarsening mechanism and indeed speeded it up so that the average radius 
eventually grew in proportion to i 1 / 2 instead of i 1 / 3 . On the other hand, if the 
precipitates were softer than the matrix, the elastic interaction was found to 
slow down the coarsening or even stop it altogether. 

A serious disadvantage of the simulations of Enomoto and Kawasaki is that 



they depend on Eshelby's formula (56) for the elastic interaction, which in turn 
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depends on an assumpton that the precipitates are spherical. As pointed out 



by Onuki and Nishimori (1991) |90| (see section 2.2.3 above) if the precipitates 
are allowed to change their shape, the interaction energy may deviate greatly 
from this formula and so these simulations do not even convincingly establish 
the i 1 / 2 growth law. 



Shortly after these simulations were done Leo et al. [ LOO | gave a general 
scaling argument leading to the conclusion that when the precipitates are large 
enough for the surface tension to be unimportant the average linear size of the 
precipitates will indeed grow in proportion to t 1 ^ 2 . Their argument is based 



on the same idea as one given by Mullins and Vihals [ 136 1 [117| which predicts 
a t 1 / 3 growth law for the case where there are no elastic forces. The main 
assumption used in this argument is that the domain structure is statistically 
self-similar, i.e. that at any time it has the same statistical properties as the 
one obtained by a uniform expansion of the domain structure at any earlier 
time. However, in contrast to the capillarity-dominated case, where there is 
plenty of evidence for the statistical self-similarity hypothesis and the t 1 ^ 3 law 
is well-established, in the case where elastic forces dominate there appear to be 
many cases, particularly for anisotropic materials, where the hypothesis of self- 
similarity is false. For example, the hypothesis would require the equilibrium 
shape of a single inclusion to be independent of its size, but for anisotropic 



elasticity this is not so (see section 2.3.1 above). The scaling argument does, 



however, tend to confirm the t 1 / 2 growth law found by Enomoto and Kawasaki 
for spherical precipitates which were harder than the matrix. 

Besides the growth law, the positional correlations of particles have received 
much attention in the non-elastic case. For elastically interacting systems, the 
main studies of positional correlations using the sharp-interface method have 
been done by Abinandanan and Johnson (1993). They calculated |33l Q the 
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concentration field in a three-dimensional isotropic matrix containing spherical 
inclusions, using a multipole expansion method. (For an alternative to the 
multipole expansion method, see ref. ]5^|). Unlike Enomoto and Kawasaki, 
they took the elastic stiffness matrix to be homogeneous, but they avoided 
the Bitter-Crum theorem by making the misfit strain tetragonal rather than 
isotropic. They obtained the effect of diffusion in the matrix (though not inside 
the precipitates) not only on the sizes of the precipitates, but also on their 
velocities. The resulting formulas were used in a simulation of the time evolution 
of a system of spherical precipitates whose initial positions and sizes were chosen 
from a random distribution. In the final configuration, the positions of most 
of the surviving precipitates were strongly correlated : arranged in a plane 
perpendicular to the axis of the tetragonal misfit strain if the axial component of 
the anisotropic part of the misfit strain was an expansion, but in lines parallel to 
the misfit strain axis if this component of the misfit strain was a compression. In 
a later paper p3[ ] Abinandanan and Johnson discuss the development of spatial 
correlations during coarsening. For small particles, capillary forces dominate 
and the particles hardly move; the spatial correlation is that each large particle 
is surrounded by a depletion zone from which it has eaten most of the solute 
material. For large particles, elastic forces dominate and the particles do move, 
setting up a new system of correlations characterized by clustering of favourably 
oriented similar sized particles. 

Hort and Johnson (1996) |^5| simulate the time evolution of a system of 
spherical precipitates in the presence of a uniaxial applied stress, both for homo- 
geneous elasticity with tetragonal misfit strain and for heterogeneous elasticity 
with isotropic misfit strain. The precipitates are found to align themselves in a 
way that is consistent with eqn ( |55| ) . 

All the calculations mentioned so far in this section depend on the (danger- 
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ous) simplifying assumption of spherical precipitates. Jou et al. (1997) |T3| do 
two-dimensional calculations which do take shape changes into account, using 



the boundary integral method mentioned in section 3.2.2. Their calculations 
refer to an elastically isotropic but heterogeneous material. Complicated effects 
are observed : alignment, movement, coagulation, coarsening, with everything 
depending strongly on the amount of elastic inhomogencity, the misfit strain, 
and the externally applied field. 

4 Diffuse interface models 

In this type of model we no longer represent the phase boundaries as a geomet- 
rical surface. Instead, the micro-structure of the interface is itself included in 
the model. The concentration c(x, t) of the A-component is now treated as a 
function of position x and time t which is continuous throughout the whole of 
VL. The surface Y and the domains f2 a and Q 13 play no part in the mathematical 
formulation of these models, though it may be possible to identify them approx- 
imately once the concentration field is known. In particular, at late stages of the 
coarsening process, large domains develop in which c is almost uniform, being 



close to one or other of the equilibrium values c" g and cf q defined in section 4.2 
these domains, which are the analogues in this model of the domains Q a and Q.P 
used in the sharp-interface model, are separated by a relatively thin transition 
layer whose thickness does not change with time, and which is the analogue of 
the geometrical surface Y in the sharp-interface model. 

An advantage of the method is that it gives information about the very early 
stages of spinodal decomposition, before any well-defined interfaces have formed 
Jl95| . For the late stages, it has further advantages: no a priori assumption 
about the topology of the interface need be made - this topology is determined 
by the evolution of the model itself - and in numerical work the difficult problem 
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of dealing with moving free boundaries is avoided. Moreover, even when the 
topology does not change, the rather complicated conditions on the surface of 
r need not be allowed for explicitly. 



4.1 The Cahn-Hilliard equation and some of its relations 
4.1.1 The Cahn-Hilliard equation 

The diffuse-interface method was introduced in Cahn's celebrated paper about 



spinodal decomposition (1961) [ 195 1 . He based his initial discussion on a for- 



mula, due originally to van der Waals [ 209 1 , for the free energy of a non-uniform 
system without elastic energy 

F = ^[/(c) + i X (Vc) 2 ]d 3 x. (78) 

Here c (meaning c(x, t)) is the concentration field, assumed to be a continous 
(indeed differentiable) function of position, /(c) has the same meaning as in (||) 
and x is a positive constant, so that the gradient term penalizes rapid spatial 
variations of c. To describe phase transitions, the function / is assumed to have 
a non-convex graph, as illustrated in Fig. 6. For example the approximation to 
/ given by the mean-field theory of a binary substitutional alloy, 

,f(c) = fi eq c - 2kT c(l - c) + kT{clogc + (1 - c) log(l - c)} (79) 

where /i eq and To are constants, has a non-convex graph when T < Tq. 
For the time dependence of c, Cahn assumed 

dC V.(W=^ (80) 



dt \ <5c(x 

where M is a mobility coefficient, which may depend on c, and 6F/5c(x.) is the 
variational derivative of F with respect to c(x). The variational derivative is 
defined by the condition that 

SF A/ ^, 3 r _ (dF{c + e(, 



x)d 3 x = v ■ (81) 
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for all differentiable functions 4> which vanish on the boundary of f2. In the 
present case this definition gives 

SF 

— =/'( c (x))- X V 2 c(x) (82) 
where /'(c) denotes the derivative of /(c), the same thing as the local chemical 



potential used in section 3J, and so flBOD becomes 
Br 

- = V-AfV(/'(c(x))- X V 2 c(x)) (83) 

This is the Cahn-Hilliard equation, which is the cornerstone of the theory of 
spinodal decomposition without elastic interactions. If the formula (|79| ) is used 
for /, M is generally taken to be Dc(l — c)/kT where D is the diffusivity 
(assumed for simplicity to be independent of c — compare eqn (|57|)), so as to 
make ( p3| ) reduce to Fick's law when c(l — c)Tq/T and x are small. 

Eqn ( |83| ) is to be solved with boundary conditions modelling the conditions 
at the surface of the material, and initial conditions modelling the initial concen- 
tration distribution. For example, to represent an alloy quenched from a high 
temperature to a point inside the miscibility gap of the phase diagram, c(x, 0) is 



generally taken (see, for example, [ 104 1 ) to be c + £(x), where c = JVa/|0| is the 
space-averaged concentration of A atoms and £ is a spatially random function 
representing small concentration fluctuations. If the uniform state c(x) — c is 
unstable (see section [4.2| below) these fluctuations grow with time, eventually 
producing concentration distributions very similar to the ones that are observed 
experimentally in spinodal decomposition if elastic interactions are unimportant. 

For many alloys, further order parameters besides the concentration are 
necessary to characterize the phase transition completely. The prototype of 
such situations is the order-disorder transition at constant concentration, as in 
Cu-Au, where at low temperatures each kind of atom prefers one of two sub- 
lattices of a body-centred cubic lattice. The mathematical model uses an order 
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parameter field ry(x) which is the local concentration of (say) Cu atoms on one 
of the sublattices. The free-energy functional has the form (|7|) with rj replacing 
c, but since r\ does not obey a conservation law the kinetic equation is not (p0[) 



but the Allen-Cahn [156 1 equation 



dv i SF 

i = - a ^r (84) 

(also known as the time-dependent Ginzburg-Landau equation) where A is a 
kinetic coefficient, which may depend on 7/. In general, both the concentration 
field and one or more order parameter fields will appear in the free-energy 
functional, the concentration obeying an equation of the form ( |S3| ) and each 
non-conserved order parameter obeying an equation of the form ( |S4| ) P| 

A further refinement of these equations, which makes it possible to allow 
for fluctuation effects such as nucleation, is to include noise terms on the right 



sides of ( |83| ) and (|S4|). This development was initiated by Cook 1 178 1 , The 
resulting equations might be called the stochastic Cahn-Hilliard and stochastic 
Allen-Cahn equations. For further information about their use in physics see the 



reviews by Gunton et al. [ 147 1 or Gunton and Droz [ 148 1 . On the mathematical 
side, it is known that the one-dimensional stochastic Allen-Cahn equation with 
space-time white noise has a solution (indeed a unique solution), but in more 
than one dimension this is known only for coloured noise. || 

4.1.2 Including elastic energy in the Cahn41illiard equation 

The free energy formula ( |78| ) can be generalized to include elastic interactions 
if we introduce a new field variable u(x) and replace / by / + w just as we did 
in (|J). In order for the elastic interaction to have any effect on the diffusion, 
however, we must allow Xij mn and/or to depend on c, so that the elastic 

11 In principle, cross-terms might appear in these equations, for example, if there are two 
order parameters, drji/dt = —An(6F/Sr]i) + Ai2(5-F/<5»?2))) 



5G 



energy density w now depends on c as well as the strain. The free energy 
formula thus becomes 

F{c(0, u(-)} = _£[/«x)) + w(c(x), eil (x), . . .) + ix(Vc(x)) 2 ]d 3 x (85) 

with u>(c, en, ei2, . . . , 633) given by the formula (|5|), as before. The c-dependence 
of w arises from the possibility that Xij mn and e^- may depend on c. Cahn's 
recipe (|o|) now gives, in place of (p^) , 

9c „ ( t i,„,„^ , 9w(c(x),en(x), . . .) 2 



- = V ■ MV (/(c(x)) + ^ ^ ^ - X V 2 c(x) 

= V • MV(/t - xV 2 c(x)) (86) 

where /2 is the diffusion potential, a function of c and en, ei2, . . . defined by 

d(/(c) + w(c,en,ei2, ■ ■ ■)) 
dc 

/i(c) + 1 ^ ^^iL (Ae) .. (Ae)mn _ £ X .. mn{Ae) ..^k± (87) 

ijmn ijmn 



Equation (|86|), without the gradient term, is due to Larche and Cahn [152|. The 



complete equation was given by Onuki (1989) [HE 



In using eqn (86) the elastic field variable u is assumed to relax to its equilib- 
rium value much faster than the concentration, and therefore it is determined, 
just as in the sharp-interface model, by minimizing F. The Euler-Lagrange 



equation for the minimization is, as in the sharp interface model, eqn (11). In 
consequence of this minimization, the strain at each point in the material is 
a (non-local) functional of the entire density profile c(-) and so the diffusion 
potential at that point is likewise such a functional. For an example, see eqn 
© below. 

As in the sharp interface model, the minimization to find u is easiest when 
the elastic stiffness tensor is independent of c and therefore uniform, so that the 



Fourier transform formula (38) can be used. A convenient way to allow for the 
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dependence of stress-free strain on concentration in this formula is to assume a 
linear relation analogous to Vegard's law : 

e% = CHj + b ljC (88) 

where a%j and bij are consants which can be expressed in terms of the quantities 
efj, c" q etc. used in chapter 2. By substituting this into eqn ( |38| ) and using the 
result in (g^) we can express the free energy as a functional of the field c alone : 

1 ^,3. 1 



F{<-)} = yj/( C (x)) + - X (Vc(x))^}d J x+- J y^( C (x)-c)^(x-y)(c(y)-c)d^d d y 

(89) 

where c is the space average of c, and V e i is defined by 

V e i(z) = ^2 hjHij mn (z)b mn (90) 

ijmn 

Apart from a numerical factor [c eq }~ 2 , V e i is the inverse Fourier transform of 
the function B in Khachaturyan's formula (|36|). The diffusion potential at x 
when F is given by (|8^) is the functional derivative of the non-gradient part of 
F with respect to c(x) , i.e. 

A(x) = /i(c(x)) + f V el (x - y)(c(y) - c) d 3 y (91) 
Jn 

wwhere /x(c) = f'(c) is the chemical potential as in ( |60| ) 

To obtain information about coarsening from such models, one must inte- 
grate (|8q) over long time intervals. Methods of integrating the non-elastic ver- 



sion of this equation are reviewed in |106|. The elastic calculation is a formidable 



one, particularly in three dimensions, since to find dw/dc the equation of elastic 
equilibrium ( pi] ) must be solved at each time step — either as a partial differ- 
ential equation in position space, or (if the elastic constants are independent of 
concentration) by the Fourier transform formula ( |3l| ) or an equivalent position- 
space formula such as the double integral in (|S9"|). Until recently it was not 
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possible to calculate time evolution from these equations with any accuracy in 
a reasonable time, and inaccurate or uncertain methods had to be used instead 
[126, 108 1; but more reliable methods are now becoming available B [?]]. 



In the limit where \ ~ * 0> with a suitable concurrent re-scaling of the time 
variable, the Cahn-Hilliard equation (|S3|) converges to the non-elastic version of 



the diffusion problem with moving interfaces, formulated in chapter 3.1. This 



was shown formally by Pego [12C] and rigorously by Alikakos et al pi] . The 
corresponding convergence problem for the elastic version, eqn (|S6|), is treated 
by Leo et al. (1998) § 

4.1.3 Lattice differential equations 

The Cahn-Hilliard equation does not take explicit account of the lattice struc- 
ture of the crystal. In the case where the elastic stiffness tensor is uniform, it 
is possible to take account of this lattice structure, without departing from the 
general philosophy of the method, by replacing the integral in ( |89| ) by a sum 
over N lattice sites situated at postions r, so that the free energy F is no longer 
a functional but a function of N variables c(r): 

F = -^E /M r )) + \ E E v ( r - r ')( c ( r ) - s ) W r ') - 5 » ( 92 ) 

r r r' 

where V is defined by 

V(a)=A(z) + V r aJ (a) (93) 

with A(z) chosen so that A(r — r')c(r') is a finite-difference approximation 
to x^ 2 c( r )- F° r example, on a simple cubic lattice with lattice spacing a, we 
could take A(r) = x/° 2 f° r the 6 lattice points nearest to the point r = 0, 
A(0) = — 6x/a 2 , and A(r) = for all other lattice points. However, once the 
free energy has been written in the discrete form (^2|) , new possibilities become 
available, since the non-elastic forces can now be allowed for explicitly instead 
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of through their contribution to the thermodynamic free energy density /(c). 
For example, V(z) could be chosen to be negative (favouring unlike pairs of 
atoms) if the vector z connects a nearest-neighbour pair of sites but positive 
(favouring like pairs) if z connects a next-nearest neighbour pair. Free-energy 



formulas of this type are the basis of Khachaturyan's discussion [194, 190, 191 1 
of equilibrium structures in crystals. 

The discrete version of the Cahn-Hilliard kinetic equation ( p0| ) is 

dc(r,t) OF dF \ 

-yj=^( r -o^___j (94) 

where F is given by (|92|) and the translationally invariant 'mobility matrix' 
M(r — r') takes a value proportional to the mobility M used in (pCj). Eqn 



( P4| ) is due to Khachaturyan [184, |150| and, independently, Hillert, Cook et 
al. [197, 177]. It was introduced into this field by Chen and Khachaturyan 

A difficulty with most derivations of (|94|), including the one just given, is that 
it is not obvious how M(r— r') should depend on c(r) and c(r') (corresponding to 
the factor c(l — c) in the continuous-space formula M = Dc(l — c) /kT given just 



after eqn (|83j)). There is a microscopic derivation (Binder (1974) [167], Penrose 
(1991) fl93[ ) leading to a more complicated kinetic equation which, if c(r) does 
not vary too rapidly from one lattice site to the next, reduces to ([)4|) with 
M(y — r') approximatetely proportional to Dc av (l — c a „)/kT where c av is the 
local average ^(c(r) + c(r')). In actual computations using (|94| ) with the elastic 
term M(r — r') has, however, been given the value Dc(l — c) /kT, independent of 
the local value of c. The results of these computations are discussed in section 



4.3 below. 



4.2 The stability of a uniform mixture 



GO 



4.2.1 Isotropic elasticity 



It was pointed out by Cahn (1961) [ 195 1 , that clastic interactions tend to stabi- 



lize a homogeneous mixture which on purely thermodynamic grounds might be 
expected to separate into two phases. In his original paper he showed that this 
mechanism would work even in an isotropic solid F^L 

The mixture is assumed to obey Vegard's law, that the stress-free strain is 
isotropic and varies linearly with the concentration, so that 

e i 3 = V( c - c o)% (95) 

where r\ and cq are constants. Then we can insert the formula (|3^) for the total 
elastic energy into (|85| ) and obtain, in the absence of externally applied strain, 

F = X (/(c) + r^; (c(x) _ 5)2 + ^( Vc ) 2 ) d3x ( 96 ) 

where c = f2 _1 J n c(x) d 3 x. The uniform state with c(x) = c everywhere is 
stable if F is a minimum with respect to variations of c satisfying the constraint 
J(c — c) d 3 x = 0. A sufficient condition for such a minimum is for the function 
/(c) + (2E/(1 — v)){c — c) 2 to be locally convex near c — c, and if 2E/(1 — v) 
is large enough this can happen even though f(c) is not locally convex, so that 
without the elastic interaction the mixture would be unstable. The uniform 
state is stable against all such variations if and only if the function /(c) + 
(2E/(1 — v))(c — c) 2 is convex, and it is stable against small variations if and 
only if this function is locally convex near c = c, that is to say if 

/"(c) + 2E/{1 -u)>Q (97) 



12 The 'Bitter-Crum' theorem might seem to ilmply that isotropic homogeneous elasticity 
never has any effect, but this is not so: what it really tells us is that in a mixture of two 



phases F depends only on the overall mole fraction, not on the geometry [ 142 ] , a fact which 
is perfectly consistent with eqn (pd)- 
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where /" denotes the second derivative of the function /. The elastic term 
can make the left-hand side of (^) positive even though /"(c) is negative, and 
so the coherent mixture can be stable against phase separation even though 
an incoherent mixture of the same composition would not be. This effect can 
lower the critical temperature of such mixtures significantly. ^] For a discussion 
of this and related effects using a simple thermodynamic model see Cahn and 
Larche (1984) p§ 



4.2.2 Anisotropic elasticity 

The above treatment of stability is easily generalized to anisotropic elasticity 



(Cahn (1962) [193|). Substituting Vegard's law (|9J) into (|3l|) we obtain, in the 



absence of externally applied stress, 



' k^O im 

n 2 

^ 2jb E l c ( k )! 2 ^E *~(k) 

' k^O im 



= -rj 2 minV^ ^ii, nm (k) [ (c(x) — c) 2 d 3 x by Parseval's theoren(i98) 

2 k tr J ° 

As in the discussion of (|96|), it follows that the uniform state is stable against 
small concentration fluctuations if 

/"(c) + -^min^* 

iimm 

(k) > (99) 



mi 



13 Khachaturyan, on page 490 of reference |150fl contests this conclusion and gives a formula 
implying that for an isotropic crystal the stabilizing effect is zero. The difference of opinion 
appears to arise from different definitions of the elastic energy: according to our definitions 
the elastic energy of an alloy without externally applied stress is zero if the concentration is 
uniform, since the mininum of J tod 3 x is zero, achieved by taking e;j = e?^. everywhere (see 
eqn (^)); whereas according to Khachaturyan's definition part of the energy is elastic energy 
even when the concentration is uniform. 
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This equation, in the form appropriate to a cubic crystal, is due to Cahn [ 193 1 . 
The elastic term is always positive (see the text just after eqn ( |33| ) above) and 
therefore the elastic interactions tend to stabilize the mixture against phase 
separation, just as they do for isotropic elasticity. 

A generalization of this kind of stability analysis, with the elastic moduli 



depending on concentration, has been used by Thompson and Voorhees 19 to 
discuss the effect of an externally applied stress on orientation of the newly- 
formed precipitates in an anisotropic alloy. 

4.3 Simulations with the diffuse interface model 
4.3.1 A single precipitate 

Wang et al. |)5| carried out two-dimensional calculations for a single precipitate 
in a cubic solid with negative isotropy, using Khachaturyan's kinetic equation 
(J)4|). The initial configuration was chosen to represent a single circular precip- 
itate, and they waited for an equilibrium to be reached, for various different 
values of a parameter representing the ratio of the elastic to the short-range 
interactions in their model. Allowing for the difference between two and three 
dimensions, this parameter is proportional to the one denoted by l/i?o m section 



1.1 of the present paper. As the parameter was increased (which corresponds to 
considering larger and larger precipitates in a real physical system) the equilib- 
rium shape, at first round, became more and more square; then, a hole appeared 
in the middle of the square and finally it split into two rectangles. It is also 
possible for the faces of the 'square' to become concave [[l8) 

In three dimensions there are even more possibilities, particularly if the 
transformation strain is anisotropic (Wang and Kachaturyan 1996 |Q). How- 
ever, as these authors show in |5(| the behaviour of isolated precipitates is not 
necessarily a reliable guide to their behaviour in the presence of others. 
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4.3.2 Many precipitates 



Nishimori and Onuki (1990) 1 101 ] carried out two-dimensional simulations us- 
ing the elastic Cahn-Hilliard equation (^6|) with an anisotropic stiffness tensor. 
Starting from an initial concentration distribution corresponding to a homoge- 
neous disordered state (i.e. a constant plus a small spatially random term), they 
found that the shapes of the precipitates were strongly affected by the elastic 
anisotropy, and that they tended to align themselves along the crystal axes. See 
Fig. 7 d,e. 

Wang et al. |)6| ^2|, |7lJ (see also ^2|), also did simulations on this type of 
system, but they used the Khachaturyan equation (^) and paid more attention 
to the effect of volume fraction and to to the visible details of the evolution 
process. For small volume fractions, the latest observed state consisted mostly 
of nearly round precipitates near the sites of a a square superlattice. This con- 
figuration was reached by a form of inverse coarsening, in which the precipitates 
that were close to sites of the nascent superlattice grew at the expense of those 
that were not, regardless of the initial sizes of these precipitates. For larger 
volume fractions (say 50 %) the precipitates, initially disconnected, tended to 
coalesce, and the latest observed state consisted of needle-like precipitates par- 
allel to the two soft directions. In a later paper |0, they consider a case where 
the transformation strain is tetragonal rather than isotropic (modelling Mg- 
stabilized cubic Z1O2): intitially, the precipitates line up along crystallographic 
directions; later, an alternating band structure develops as the precipitates in 
some rows grow at the expense of those in others. 



Onuki and Nishimori (1991,1992) |92j |7§] were the first to investigate the 
effects of elastic inhomogeneity by means of this type of simulation. To isolate 
the effect of the inhomogeneity (which, confusingly, they call 'elastic misfit' in 
some of their papers) they assumed isotropic elasticity. If the softer (less rigid) 
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phase occupied less than half of the total volume, then at early times the config- 
uration consisted of isolated precipitates of softer phase embedded in a matrix 
of the harder phase, but as coarsening proceeded these precipitates joined up to 
make a percolated network structure, so that it was now the harder phase which 
was arranged in isolated inclusions, while the softer phase formed the surround- 
ing matrix (see Fig. 7(a-c)). Moreover, the growth of these precipitates was 
anomalously slow, i.e. slower than the usual i 1 / 3 law. Later, these authors did 
simulations of the same type for anisotropic (but still inhomogeneous) elasticity 
[ p9| . The anisotropy made the inclusions align themselves parallel to one of the 
crystal axes, but still the softer phase wrapped the harder one. If an external 
stress in an oblique direction was applied, the preferred directions of the inclu- 
sions tilted so as to be more nearly perpendicular to one of the principal axes 
of the external stress. 

The numerical method used by Onuki and Nishimori for the time evolution 
was not an accurate one; later work using accurate methods (Leo et al. 1998 [Q) 
has, however, confirmed their conclusions about the reversal of the topological 
roles of the two phases. 

These topological changes, which the sharp-interface model had completely 
failed to predict, can be understood qualitatively on the basis of the formula for 
spherical inclusions given in section ( [2.2.1 ) , according to which the energy per 
unit volume of precipitate is less for rigid precipitates in a flexible matrix than 
for flexible precipitates in a rigid matrix. Or, even more qualitatively, consider 
the extreme case where both phases have the same bulk modulus but one of the 
phases is perfectly flexible, i. e. its shear modulus is zero. If the inclusions are 
perfectly flexible the rigidity of the surrounding matrix will prevent them from 
taking up the sizes where their density has its stress-free value and therefore 
they will have some strain energy; but if the matrix is perfectly flexible then the 
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(rigid) inclusions can take up their stress-free sizes and shapes and the matrix 
can expand or contract to its stress-free density, so that there is no strain energy 
at all. Thus the energy is lower for a perfectly flexible matrix than it is for 
perfecty flexible inclusions. 

4.3.3 Ordered phases 

In alloys such as Ni-Al one of the phases (in this case the M3-AI phase) can 
have more than one variant. Such systems can be treated using Khachaturyan's 
kinetic equation ( |94| ) , generalized by incorporating additional order parameters 
as described at the end of section |4.1.l| . Simulations of various types of alloy 
requiring this type of treatment have been done by Sagui et al (1994) [|50| , by 
Wang and Khachaturyan (1995) fc} , by Li and Chen (1998) §,[l(J and by Wang 
etal. (1998) Jl§]. An important new possibility is that the transformation strain 
can be anisotropic and hence different in the different variants; then, when an 
external stress is applied, one variant will grow at the expense of another. In 



this way, since the conditions of Eshelby's theorem (see section 2.1.4) do not 
apply, rafting can take place even with homogeneous elasticity. 

A recent paper of this type by Le Bouar et al. |l[] models the time evolution 
of alloys in which the matrix is a cubic crystal but the precipitates are tetrago- 
nal, using a generalized stochastic Cahn-Hilliard equation with three dependent 
variables : the concentration and two order parameters. The theoretical model 
is spectacularly successful : see Fig. 2. 

4.3.4 Thin films 

In this section we look at the possibility of elastic interactions arising from a 
violation of condition 4 of the Bitter-Crum theorem, that is to say in a system 
with boundaries. An effect of this type was noted by Larche and Cahn |[5| 
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who pointed out that if spinodal decomposition takes place in in an isotropic 
crystal with uniform elastic constants the equilibrium state can have the two 
phases on opposite faces of the film, and because of the elastic misfit the film 
will bend. This effect has indeed been observed, in the a — a' phase separation 



of Nb-H 16C, 155]. In these experiments, the hydrogen-rich a' phase was found 
on the outer side of simply curved specimen foils. In some cases the foil was 
even observed to split into two pieces with half thickness. Both halves were 
strongly curved. These results were interpreted as a spinodal decomposition 
process where the specimen thickness was equal to half the spinodal wavelength 
for the simply curved foil and exacctly to this wavelength for the split foil 
Putting an appropriate approximation for the elastic strain field into the 



elastic Cahn-Hilliard equation (86) Larche and Cahn obtained the kinetic equa- 
tion 

-df(c) , 2 V 2 E ( _ &r \ ^ 



— = V 
dt 



dc 1 — v \ 2h 



(100) 



i r h 

M? = — / (c- c)xdx (101) 
« J-h 



where r] is the compositional expansion coefficient, defined in (95), E and v are 
Young's modulus and Poisson's ratio, c is the average concentration, and Mf 
is the first moment of the concentration distribution in the direction at right 
angles to the film, defined by 

-h 
i-h 

Here h is half the thickness of the film and the coordinate system in the material 
is chosen to make the two faces of the unstrained film the planes x = h and 
x = —h. The one-dimensional version of these equations, obtained by neglecting 
the dependence of c on y and z, was studied both analytically and numerically by 
Larche and Cahn, and numerically in more detail by Cahn and Kobayashi p5| . 
The most striking result was that the elastic term accelerated the separation 
into two phases. 
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Another bounded-system problem to which the diffuse-interface model can 
be applied is the growth of thin films by deposition, an important process in the 
semiconductor industry. The elastic problem is more complicated than before 
because there is now a misfit between the film and the substrate as well as the 
one due to compositional nonuniformity. Moreover, the possibility of surface 



diffusion [203, 121 and the effect of material arriving at the surface have to be 
allowed for. At the time of writing the implications of linear stability theory are 
still being worked out |(39|, |24|, ||, ||] , though nonlinear effects are clearly 
important too |7Q|. 



5 Atomic lattice models 
5.1 The Born model 

In this approach continuum elastic theory is abandoned; instead the elastic prop- 
erties of the material are obtained by considering the mechanical forces between 
individual atoms. Each atom is assiged a definite site p on the lattice, whose 
position when the lattice is undistorted we denote by x p , and the displacements 
Up of the atoms from these sites are assumed to be small enough to justify using 
the harmonic approximation, i.e. treating the force on each atom as a linear 
function of all the atomic displacements of all the atoms. The way to calculate 
the macroscopic elastic properties of a crystal from the interatomic force param- 
eters is explained in the book by Born and Huang (1954) ||204f , and, by doing the 
calculation in reverse, the interatomic force parameters can be obtained from 
phonon dispersion curves. For simplicity, we assume that only two-body forces 
act, so that atoms which are close enough together to interact can be thought 
of as being joined by springs. However, one should not be misled by this picture 
into assuming that all the atomic interactions are necessarily central forces; this 
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assumption would imply certain relations between the elastic constants (some- 
times called Cauchy relations) ]204| which are violated in many real materials. 
For example in a crystal with cubic symmetry the Cauchy relation is C\2 = C44 
which, although a good approximation for some ionic crystals, is a very bad one 
in most metals. 

In the case of an alloy, the interaction between two atoms depends not 
only on their separation but also on their nature. It is usual to think of their 
interaction energy as the sum of two parts, one of which, the so-called chemical 
interaction, depends only on what kind of atom they are and not on their 
positions; the other part, the elastic interaction, does depend on their positions 
and has a minimum value of zero which is achieved for some relative position. 
This energy-minimizing relative position depends on the nature of the two atoms 
: one can think of the two kinds of atoms as being of different sizes, and the 
natural length of the spring joining the centres of two atoms the sum of the 
radii of the atoms. The stiffness of the spring may also depend on the nature 
of the two atoms. The theory of the equilibrium properties of solid solutions 



based on this type of model has been developed by Khachaturyan [184, 150] and 



by Cook and de Fontaine (1969,1971,1972) |179j, [175|, |17J. For example Seitz 



and de Fontaine (1978) [ 163 1 have used the model to do energy calculations 
for various arrangements of a pair of precipitates in Al-Cu alloys; the results 
indicated a tendency for the precipitates to stack as parallel disks, providing 
a possible explanation of the oberved behaviour of supersaturated Al-Cu solid 
solutions. 

Diffusion can be represented in this model by allowing atoms of the two 
different kinds to exchange places from time to time, though sufficiently rarely 
to allow the lattice to come to equilibrium under the new forces after each ex- 
change. These exchanges can be simulated by the Monte Carlo method known 
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as Kawasaki dynamics, in which the following sequence of steps is carried out 
a large number of times : (i) choose a pair of neighbouring atoms at random, 
and also a random number £ uniformly distributed between and 1, (ii) cal- 
culate the increase in free energy AF that would occur if the two atoms were 
exchanged, and (iii) carry out the exchange if and only if £ < p(AF), where 
the transition probability function p(-) satisfies the detailed balance condition 
p(x)/p(—x) = exp(— x/kT). The change in free energy consists of two parts : 
one is the change in (mean) energy which is equal to the change in the minimum 
potential energy as a function of the atomic displacements, since the vibrational 
energy is an unchanging kT per degree of freedom. The other part is an entropy 
term depending on the vibrational spectrum of the lattice, which is generally 
neglected. Indeed, this term is zero if the atomic stiffness parameters are in- 
dependent of the way the two kinds of atoms are arranged on the lattice sites, 
a situation which makes the material elastically homogeneous. Just what part 
the entropy term should play in models capable of representing heterogeneous 
elasticity is unclear. 

We shall describe two types of simulation in which Kawasaki dynamics has 
been used to follow the diffusive motion of the atoms; one using central forces 
and heterogeneous elasticity, and the other using noncentral forces and homo- 
geneous elasticity. 



5.2 Central forces 
5.2.1 The model 

In two dimensions, the elastic stiffness tensor of any crystal with hexagonal 
symmetry is isotropic; moreover, if only central forces act, then there is a Cauchy 
relation connecting the elastic constants which is equivalent to v (i.e. Poisson's 



ratio) = 1/4 [169|. In three dimensions, the elastic stiffness tensor of a harmonic 
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crystal with fee symmetry and only central forces is likewise isotropic Jig ] and 
again Poisson's ratio is 1/4. Lee |4|, took advantage of these circumstances 
to simulate inclusions in a isotropic matrix using an atomic lattice model. He 
took the elastic energy to be the sum of the energies of a set of longitudinal 
harmonic springs, one for each pair of nearest neighbour sites, and both the 
stiffness and the natural length of each spring depended on the nature of the 
atoms at its two ends. In symbols, the elastic energy is 

W = \ E E L & - P'' Tp. 7p')(W - *(P - P'» 7 P , 7p')) 2 (102) 
p p' 

where sums go over all the sites on some lattice, 7 P is defined to to be +1 
if site p is occupied by an A atom and —1 if it is occupied by a B atom, 
L(p — p',7,7') and Z(p — p',7,7) are the stiffness and natural length of the 
spring connecting an atom of type 7 at site p to an atom of type 7' at site p', 
and r p p> — |x p — x p > | + (u p > — u p ) ■ (x p — x p / ) is the linearized distance between 
the atoms at sites p and p'. 

5.2.2 The results 

For a system of two inclusions in two dimensions, Lee found p2[ that if the 
inclusions were softer than the matrix they deformed and eventually coalesced; 
if they were harder, they remained circular and the larger one grew while the 
smaller one shrank to zero. This work was extended to anisotropic elasticity in 
ref p(|, with similar results. 

In ref Jla] Lee (1997) studies the effect of elastic inhomogeneity on a system 
of two inclusions for which the misfit strain and the elastic constants can both be 
different. If one inclusion is hard (i.e. harder than the matrix) and has a misfit 
strain while the other is soft one and has no misfit, or a misfit of the opposite 
sign, then at equilibrium (in two dimensions) the soft inclusion surrounds the 
hard one. But if the soft inclusion has a misfit strain and the hard inclusion has 



71 



none, the inclusions move away from one another; and if the two inclusions have 
misfit strains of the same sign, the soft inclusion partially wets the hard one 
but the rest of the soft inclusion tries to get away. Similar results are obtained 
in three dimensions. If an external stress is applied then the inclusions change 
their shape in a way compatible with eqn (|55|). Reviews of the work described 
in this section are given in and 0. 

5.3 Non-central forces 
5.3.1 The model 

Fratzl and Penrose (1994) |3?J describe a model for use in atomic lattice simula- 
tions which is intended to represent coarsening in a cubic crystal with realistic 
elastic constants. For reasons explained earlier in this chapter, it is essential to 
include some non-central forces in the model. This can done by adding to the 
right side of the 'longitudinal' elastic energy formula (|l02| ) a term representing 
transverse springs 

\ E E T (p - p'H( u p ^ u p') A ( x p - x p')/I x p - x p'|} 2 (los) 

p p' 

In these simulations, a square lattice was used, with elastic forces which can be 
pictured by imagining the atoms connected by three different types of spring : 
longitudinal springs connecting nearest neighbors, transverse springs connecting 
nearest neighbors, and longitudinal springs connecting next-nearest neighbours. 
To model the elastic misfit, the natural length of a longitudinal spring connect- 
ing sites p and p' was made to depend linearly on 7 P + 7 p / , but the stiffnesses 
of the springs were taken to be independent of what type of atoms they con- 
nected. This model is characterized by three stiffness constants, the same as the 
number of elastic constants in a cubic crystal. It is possible, in fact, to relate 
the microscopic stiffnesses to the macroscopic elastic constants. Formulae for 
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different lattices are given in, for example, 17£, 15C, |37j]. 

Given any configuration 7(p) of A and B atoms on the lattice, the minimum 
of the elastic energy with respect to the atomic displacements u can be found 
as a function of the 7(p)'s. The elastic free energy differs from this minimum 
by a constant, which can be ignored. Since the elastic constants are taken to be 



homogeneous, this calculation can be done using Fourier transforms [179, 15C|, 
in a similar way to the derivation of Khachaturyan's formula ( |36| ) given in section 
|2.1.3| (eqns @-(||)). The result can be written 

W= ^E B ( k ) l7(k)| 2 +const. (104) 

where the sum extends over the first Brillouin zone of the lattice and 7 is the 
Fourier transform of 7(p) : 

7 (k)=^7(p)e*P- k (105) 
p 

The 'elastic potential' B(k) can be calculated for any type of lattice, and has 



been given explicitly for some, e.g. in [17£, |37| . 

The elastic energy ( |l04j ) can be written in an alternative form, as a sum of 
interactions between pairs of lattice sites analogous to the double sum in eqn 
( |92| ) . Adding to it the chemical energy, the important part of which can also be 
expressed as such a sum, we obtain for the total free energy 

F = W + l -Yl (P - P') + V -hem(p - P')) 7(p) T(P') (106) 

p p' 

where Wo is a constant, V e ;(p) is the inverse Fourier transform of B(k), (in 
which we take B(0) = 0, so that J2 P ^ei(p) = 0) and V c hem(p — p') represents 
the short-range chemical interaction between two like atoms on sites p and p'. 

The effects of external stress can also be treated to some extent within 
this framework 41 The treatment requires, however, a modification of 
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the model since according to Eshelby's theorem (see section 2.1.4) no effect 
will be seen for homogeneous elasticity. This modification is to assume that the 
stiffnesses of the springs depend weakly on what kind of atom they connect, and 
to treat the nonuniformity as a perturbation — for details see ]l4] . In analogy 



with the analysis for an clastic medium summarized at the end of section 2.1.3, 
the elastic energy can be expressed as the sum of three parts, one of which is 
the non-quadratic 'Eshelby' interaction. To keep the problem mathematically 
simple, this non-quadratic interaction was neglected in jl4|] . 

5.3.2 The results 

The model just described can be used for simulations in various ways. 

When the chemical interaction is taken to be a nearest-neighbor attraction 
and the elasticity has cubic anisotropy, stripe-like domains develop under the 
influence of the lattice misfit ^3|. The stronger the misfit or the stronger the 
anisotropy of the elastic constants, the faster the anisotropy of the precipitates 
increases. The mean domain size, however, (defined as the ratio of area to 
perimeter) always increases in proportion to t 1 ^ 3 as in the case without elastic 
misfits §§. 

When uniaxial external stress is included in the model, wavy parallel stripes 
develop with an orientation depending on the direction of applied stress | fli"[ [Lff . 
Thus, the cubic symmetry in the domain morphology is broken, in agreement 



with the experiments on rafting described in section 1.2.3 Surprisingly, the mean 
domain size still increased in proportion to t 1 ^ 3 in the model. 

When the chemical interaction is taken to be a nearest neighbour repulsion 
instead of an attraction fTi"|| , precipitates with atomic order (of antiferromag- 
netic type) are formed. The qualitatively new feature compared to the previous 
cases is that the precipitates can be ordered on two different atomic sublat- 
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tices. Where two domains ordered on different sublattices join, there is an 
anti-phase boundary (APB). Such APBs have the property of being wetted by 
the disordered matrix, thereby creating new boundaries between precipitates. 
The presence of wetted APBs also influences the rate of coarsening. The re- 
sulting patterns, illustrated in Fig. 8, are extremely similar to those found for 



nickel-base superalloys, as described in section 1.2.1 



6 Conclusion 

We can now consider to what extent the theoretical models we have been de- 



scribing explain the various experimental phenomena mentioned in section 1.2 
The shapes of the precipitates 

The plate-like shape of the individual precipitates in alloys such as Al-Cu 



and Cu-Be, mentioned in section 1.2.1 is consistent with the prediction of the 



sharp- interface model (section [2. 3.1 ) that when the interface energy is negligible 
and the material is elastically anisotropic the equilibrium shapes of the precipi- 
tates are plates perpendicular to an elastically soft direction. Two-dimensional 
multi-precipitate simulations with anisotropic elasticity, using either the diffuse- 



interface model (see section 4.3.2 and Fig. 7(d,e)) or the atomic lattice model 
[ p3| , also give long thin precipitates perpendicular to an elastically soft direction. 
The dependence of the shape of precipitates on their size, (more precisely, on 



the parameter R/Rq defined in section 1.1) was investigated in the experiments 



with Ni-Al-Mo alloys illustrated in Fig. 1. The individual precipitates looked 
round when this parameter was very small, but roughly square when it was 
large. All three of the models we have described explain this dependence as 
an effect of the elastic anisotropy; it is seen particularly clearly in calculations 
using the sharp interface and diffuse interface models, where the behaviour of 
a single precipitate as a function of its size (that is, of the ratio of its elastic to 
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its interfacial energy) has been considered (see sections 2.3.1 , and fOP] ). 

The splitting of large cube-like precipitates in some nickel-base alloys (see 
Fig. 4) has also been investigated using various theoretical approaches, but the 
results are inconsistent. Sharp-interface energy calculations assuming inclusions 



of fixed shape do give this effect (sections 2.3.1 and 2.3.2 , as do diffuse-interface 



calculations of the equilibrium shape (sections 4.3.1). On the other hand the 



sharp-interface evolution calculations that have been done so far (section 3.2.2) 
do not exhibit the effect. 

The arrangement of the precipitates 

Experiments show that when the (anisotropic) elastic interaction is strong 
enough to affect the shapes of the precipitates, these precipitates tend to align 
along certain directions (see Fig 1(c)). A similar type of ordering is found in 
the sharp-interface model (sections |2.3.3 , 3.2.3 , 3.2.4 ) and the diffuse-interface 
model (section |4~3^) . 



The rate of coarsening 

While the t 1 / 3 growth law for coarsening in alloys without elastic misfit is 



reasonably well understood [198, 199, 117], theory has been less successful in ex- 
plaining the growth law when elastic forces do act. The experiments illustrated 
in Fig. 3 indicate [Q that if the elasticity is strongly anisotropic and (rela- 
tively) weakly heterogeneous, the coarsening follows a t 1 ' 3 growth law, whereas 
if it is weakly anisotropic and strongly heterogeneous, the coarsening follows a 
t 1 ! 3 law only for a limited time, and then slows down or stops. The atomic 
lattice simulations described in Section 5.3.2 also give t 1 ! 3 growth for a purely 



anisotropic lattice, but theory does not explain why the elastic forces apparently 
do not affect the growth exponent at all in this case. For a heterogeneous lattice 



without anisotropy the diffuse-interface simulations described in Section 4.3.2 
predicted that the growth would be slower than t 1 / 3 but not that the growth 
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exponent would change rather abruptly as shown in Fig. 3(b,c). Theory also 
predicts an interesting topological rearrangement if in the initial configuration 
the precipitates are softer (less rigid) than the matrix, but these experiments did 
not reveal any striking difference of behaviour between precipitates that were 
much softer than the matrix and ones that were harder. 
Rafting 

The effect of an externally applied stress on the shapes and arrangement of 



the precipitates, described in section 1.2.3 and illustrated in Fig. 1(g), can be 
understood in principle on the basis of isotropic (but inhomogeneous) elasticity 



using the sharp interface model, as discussed in section 2.2.1. The effect has 
been studied in a more detailed way, both for isotropic and anisotropic elasticity, 



using all three of the theoretical approaches we have described (sections 2.2.2, 



4.3.2 , 5.2.2 , 5.3.2 and Fig. 8(c,d)). All these approaches reproduce, to a greater 
or less extent, the qualitative experimental features of rafting, summarized in 
equation (|55|). 

Generalizations of the model 

Figs 2 and 8 illustrate two cases that go beyond the simple binary alloy with 
either isotropic or homogeneous anisotropic elasticity which underlies most of 
the theoretical work described in this article. In Fig. 2 the precipitates have a 
different crystal symmetry from the matrix, and in Fig. 8 the precipitates can 
occur in more than one variant. In both cases, the results illustrate the excellent 
qualitative agreement with experiment that simulations based on a well-chosen 
theoretical model can give. 
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8 List of symbols 



Symbol 

B(k) 

B(k) 

Cn, C12, C44 

c 

c 
CO 

D 

E 

E(p) 



first appearance meaning 



F 
F r 

/ 
G 

z 

K 

k 
I 
i 

M 
n 
P 
P 

7? 

r p,p' 
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sec. 3.1 
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i), © 



(ED 
©H 

various 

© 

section 5.1 



section [O] 
section [Hi] 

PD 



'clastic potential' 

elastic moduli of cubic crystal 

concentration of A (solute) atoms 

zero-temperature equilbrium value of c in Q a 

concentration at which stress-free strain is 

diffusivity 

Young's modulus 

elastic two-site interaction 

strain tensor 

stress-free strain tensor 

'externally applied strain' 

stress-free strain in phase a 

free energy of system 

free energy of interface 

thermodynamic free energy density 

shear modulus 



bulk modulus 

Boltzmann's constant 

stiffness of a longitudinal spring 

natureal length of a longitudinal spring 

mobility coefnent or matrix 

a unit vector 

ex(g@s pressure inside inclusion 
a lattice site 

isotropic stress-free strain 
radius of (spherical) precipitate 
'crossover' value of R 

linearized distance between atoms on sites p, p' 



u 


(H) 

VI 1/ 


displacement field 


u' 


(PI) 


periodic part of u 


V 


© 


total effective two-site interaction 




v h 1/ 


chemical interaction 


Vet 




two-point elastic interaction 


v n 


PI) 


normal velocity of interface 


W 




elastic energy of the specimen 


w ext 


(1) 


potential energy of loading mechanism 


Wint 


(PI) 


elastic interaction energy 


w 


(PI) 


elastic energy density 


w a 


(11 


elastic energy per unit volume of inclusion 


X 


(PI) 


reference position of material point 




(PD 


Cartesian components of x (i = 1, 2, 3) 


a,/3 




the A-rich, B-rich phases 


r 


(i 

l_T 


the interface between phases a and (3 




(Col) 


+1(— 1) for site occupied by an A(B) atom 


A(z) 


(p3|) 


X times finite-difference Laplacian 


(Ae)a 


(1) 

V LJ- / 




8 A A 


(0) 


1 if i = j, if not 


7} 


PD 


compositional expansion coefficient 


V fx) 




order-parameter held 


0(k) 


(1) 






(p6h 


mean curvature of T 


A 


(PI 


kinetic coefficient in Allen-Calm eqn 




© 


elastic (stiffness) tensor 


// 


dsob 


chemical potential 


A 




dittusion potential 


// 


© 


Poisson's ratio 


7T 


(ID 


grand canonical pressure 




(D 


interfacial energy (surface tension) 


X 


© 


coefficient of gradient term in free energy density 


^ ijmn 


© 






(1) 


region occupied by specimen, inclusions, matrix 


50 


(!) 


surface of Q 



9 Figure captions 

Fig. 1: (a) Transmission electron micrograph (TEM) of Ni-Al-Mo alloy with 
Mo-composition chosen such as to make the lattice spacing in matrix and pre- 
cipitates (gamma-prime phase) equal. Treatment: 430h at 1048 K after quench 
from the single-phase region. The plane of observation corresponds to the crys- 
tallographic plane (001). One observes round precipitates i.e. no elastic effects. 

(b) Small-angle X-ray scattering (SAXS) data from the same specimen in the 
same orientation. One observes an isotropic scattering pattern. 

(c) TEM of the (OOl)-plane of a Ni-Al-Mo alloy with Mo-composition chosen 
such as to make the lattice spacing in the precipitates (gamma-prime phase) 
larger by 0.4precipitates, aligned along the elastically soft directions, [010] and 
[100]. 

(d) SAXS from the specimen in (c). There is a flower-like pattern with four- fold 
symmetry, the elongations being in [010] and [100]-directions (indicating flat 
interfaces and strong correlations in those directions). 

(e) TEM of Ni-Al-Mo alloy with Mo-composition chosen such as to make the lat- 
tice spacing in the precipitates (gamma-prime phase) smaller (misfit -0.5qualitatively 
similar to those in (c). 

(f) SAXS from the specimen in (e). 

(g) Same alloy and heat-treatment as in (e), but now with a external compres- 
sive load of 130MPa applied to it along the vertical [010]-direction. 

(h) SAXS from the specimen in (g). Note that the horizontal streak has disap- 
peared, corresponding to the disappearence of vertically oriented interfaces in 

(g)- 

The data in (a)-(d) are taken from |3(| and in (e)-(h) from Jl4] . 

Fig. 2: (a) Chessboard- like microstructure observed by transmission electron 
microscopy in Co39.5-Ptgo.5 slowly cooled from 1023K and aged 15 days at 873K. 
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The scale bar represents 30 nm. The cubic phase appears white while the 
tetragonal phase is black, (b) Computer simulation using a diffuse interface 
model (From Q). 

Fig. 3: Time-dependence of the mean precipitate size for (a) Ni-Al-Mo alloys 



(from |6|), (b) Ni-Cu-Si alloys (from |l3j) and (c) Ti-Mo alloys (from |5j). 
The coefficients M(c) in (a) and M(T) in (c) are used to scale data obtained for 
different alloy compositions or different temperatures, respectively. 

Fig. 4: Typical example of ordered precipitates in a nickel-base superally, 
which split into several pieces during the coarsening process (from [^9i ) . 



Fig. 5. Geometrical objects used in section 2.1.2 



Fig. 6. A typical free-energy density function /. 

Fig. 7: Typical snapshot pictures for isotropic but heterogeneous elasticity 
(a-c) and homogeneous elasticity with cubic anisotropy (d-e) using equations of 
the type (^5| - |8rj ) . In (a-c) , the stiffer phase is shown white with volume fractions 
of 0.3, 0.5 and 0.7, respectively (from |p2[). Note that the softer phase always 
wraps stiffer particles. In (d-e), the volume fraction of the white phase is 0.5 
and 0.7, respectively (from [ |101| ). 

Fig. 8: (a) and (c) are reproductions of Figsl(e) and 1(g), respectively, 
(b) shows results from computer simulations of an Ising model with elastic inter- 
actions on a simple square lattice and with repulsive interaction of like atoms 
on nearest neighbor sites (with interaction energy J) and attractive interac- 
tion of like atoms on next nearest neighbor sites (energy J/2). Black shows 
the disordered phase (containing mostly A-atoms) and white the ordered phase 
(consisting of about half A and half B atoms). The overall concentration of 
B-atoms was 0.35 and the run was performed at a temperature of T=0.567J/k 
on a 128 x 128 lattice and stopped after 10 6 MCS . 

(d) The same model, temperature and annealing time as in (b), but with an 
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additional external load along the vertical direction (simulations from |Tl| ) 
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